path: root/common/utils.c

#include <stdlib.h>
#include <math.h>

#include "utils.h"

/* http://www.evanmiller.org/mathematical-hacker.html */
static long int _fib(unsigned long int n)
{
    return lround((pow(0.5 + 0.5 * sqrt(5.0), n) -
                   pow(0.5 - 0.5 * sqrt(5.0), n)) / sqrt(5.0));
}

/* http://blog.noblemail.ca/2013/01/on-calculating-fibonacci-numbers-in-c.html */
long int fib(unsigned long int n)
{
    return lround((pow(0.5 + 0.5 * sqrt(5.0), n)) / sqrt(5.0));
}

/* factorial */
long int fac(unsigned long int n) {
    return lround(exp(lgamma(n+1)));
}

list_t *list_append(list_t *list, void *value)
{
    list_t *ret, *curr;

    ret = malloc(sizeof(list_t));
    ret->value = value;
    ret->next = NULL;

    if (!list)
    {
        return ret;
    }

    curr = list;
    while (curr->next)
    {
        curr = curr->next;
    }

    curr->next = ret;

    return list;
}

list_t *list_append_int(list_t *list, int i)
{
    int *p = malloc(sizeof(int));
    *p = i;
    return list_append(list, p);
}

list_t *list_append_long(list_t *list, long l)
{
    long *p = malloc(sizeof(long));
    *p = l;
    return list_append(list, p);
}

int list_len(list_t *list)
{
    int i = 0;

    while (list)
    {
        list = list->next;
        ++i;
    }

    return i;
}

void *list_get_n(list_t *list, int n)
{
    int i = 0;

    for (; list && i < n; ++i)
    {
        list = list->next;
    }

    if (list)
        return list->value;

    return NULL;
}

void list_append_list(list_t *list, list_t *elem)
{
    while (list->next)
    {
        list = list->next;
    }

    list->next = elem;
}

void list_sort(list_t **list, list_compare_cb cmp)
{
    int len = list_len(*list);
    list_t *pivot = *list;
    list_t *less = NULL, *greater = NULL;
    list_t *curr;

    if (len <= 1)
        return;

    curr = pivot->next;
    pivot->next = NULL;

    while (curr)
    {
        list_t *next = curr->next;
        curr->next = NULL;

        if (cmp(curr->value, pivot->value) <= 0)
        {
            if (!less)
                less = curr;
            else
                list_append_list(less, curr);
        }
        else
        {
            if (!greater)
                greater = curr;
            else
                list_append_list(greater, curr);
        }

        curr = next;
    }

    list_sort(&less, cmp);
    list_sort(&greater, cmp);

    list_append_list(less, pivot);
    list_append_list(less, greater);

    *list = less;
}

void list_free_int(void *value)
{
    free((int*)value);
}

void list_free_long(void *value)
{
    free((long*)value);
}

void list_free(list_t *list, list_free_cb cb)
{
    list_t *tmp;

    if (!list)
    {
        return;
    }

    while (list->next)
    {
        cb(list->value);
        tmp = list;
        list = list->next;
        free(tmp);
    }

    cb(list->value);
    free(list);
}

list_t *prime_factors_naive(const long *prime_table, long prime_table_size,
    long number)
{
    long prime = 0;
    long i;
    char quot_is_prime;

    list_t *factors = NULL;

    // printf("table size: %ld\n", prime_table_size);

    while (prime < prime_table_size)
    // while (1)
    {
        // printf("prime: %ld\n", prime);

        long quot = number / prime_table[prime];
        long rem = number % prime_table[prime];

        if (rem == 0)
        {
            // printf("factor: %ld\n", prime_table[prime]);

            factors = list_append_long(factors, prime_table[prime]);

            if (quot == 1)
            {
                break;
            }

            quot_is_prime = 1;

            for (i = 2; i < sqrt(quot); ++i)
            {
                if (quot % i == 0)
                {
                    quot_is_prime = 0;
                    break;
                }
            }

            if (quot_is_prime)
            {
                // printf("factor: %ld\n", quot);
                factors = list_append_long(factors, quot);
                break;
            }

            number = quot;
            prime = 0;
        }
        else
        {
            prime++;
        }
    }

    return factors;
}

long int gcd(long int a, long int b)
{
    if (a == 0)
        return b;

    if (b == 0)
        return a;

    while (1)
    {
        long int q = a / b;
        long int r = a % b;

        if (r == 0)
        {
            return b;
        }

        a = b;
        b = r;
    }
}

long int lcm(long int a, long int b)
{
    a = a / gcd(a, b);
    return a * b;
}

int is_coprime(long int a, long int b)
{
    if (gcd(a, b) == 1)
    {
        return 1;
    }

    return 0;
}

int pythagorean_triplet(int p, int q, int *a, int *b, int *c)
{
    // http://www.friesian.com/pythag.htm

    // Two cases:
    //
    // 1. p and q must both be odd, p>q, gcd(p,q)=1
    // 2. p and q must not both be odd, p>q, gcd(p,q)=1

    if (p > q && gcd(p, q) == 1)
    {
        if (p % 2 != 0 && q % 2 != 0)
        {
            *a = p*q;
            *b = (p*p - q*q) / 2;
            *c = (p*p + q*q) / 2;

            // sometimes a > b, we can just switch
            if (*a > *b)
            {
                p = *a;
                *a = *b;
                *b = p;
            }

            return 1;
        }
        else if ((p % 2 == 0 && q % 2 != 0) || (p % 2 != 0 && q % 2 == 0))
        {
            *a = 2*p*q;
            *b = p*p - q*q;
            *c = p*p + q*q;

            // sometimes a > b, we can just switch
            if (*a > *b)
            {
                p = *a;
                *a = *b;
                *b = p;
            }

            return 1;
        }
    }

    return 0;
}