/* 
 * @brief    Cyclical figurate numbers
 * 
 * Project Euler (projecteuler.net)
 * 
 * Problem 61:  Triangle, square, pentagonal, hexagonal, heptagonal, and
 *  octagonal numbers are all figurate (polygonal) numbers and are generated
 *  by the following formulae:
 *  
 *  Triangle        P(3,n)=n(n+1)/2       1, 3, 6, 10, 15, …
 *  Square          P(4,n)=n²             1, 4, 9, 16, 25, …
 *  Pentagonal      P(5,n)=n(3n−1)/2      1, 5, 12, 22, 35, …
 *  Hexagonal       P(6,n)=n(2n−1)        1, 6, 15, 28, 45, …
 *  Heptagonal      P(7,n)=n(5n−3)/2      1, 7, 18, 34, 55, …
 *  Octagonal       P(8,n)=n(3n−2)        1, 8, 21, 40, 65, …
 *  
 *  The ordered set of three 4-digit numbers: 8128, 2882, 8281, has three
 *  interesting properties.
 *  
 *  1. The set is cyclic, in that the last two digits of each number is the
 *     first two digits of the next number (including the last number with
 *     the first).
 *  2. Each polygonal type: triangle (P(3,127)=8128), square (P(4,91)=8281),
 *     and pentagonal (P(5,44)=2882), is represented by a different number in
 *     the set.
 *  3. This is the only set of 4-digit numbers with this property.
 *  
 *  Find the sum of the only ordered set of six cyclic 4-digit numbers for
 *  which each polygonal type: triangle, square, pentagonal, hexagonal,
 *  heptagonal, and octagonal, is represented by a different number in the set.
 *  
 */
#include <assert.h>
#include <stdbool.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

struct Node;

typedef struct {
    uint16_t length;
    uint16_t capacity;
    struct Node *items[];
} NodeArray;

typedef struct Node
{
    struct Node *next;
    uint8_t side_count;
    uint16_t number;
    uint8_t incoming, outgoing;
    NodeArray *neighbours;
} Node;

typedef uint16_t (*Function)(uint16_t n);

uint16_t triangle(uint16_t n) { return n * (n + 1) / 2; }
uint16_t square(uint16_t n) { return n * n; }
uint16_t pentagonal(uint16_t n) { return n * (3 * n - 1) / 2; }
uint16_t hexagonal(uint16_t n) { return n * (2 * n - 1); }
uint16_t heptagonal(uint16_t n) { return n * (5 * n - 3) / 2; }
uint16_t octagonal(uint16_t n) { return n * (3 * n - 2); }

static const Function FUNCTIONS[] = {
    triangle, square, pentagonal, hexagonal, heptagonal, octagonal
};
#define FUNCTION_COUNT  (sizeof(FUNCTIONS) / sizeof(FUNCTIONS[0]))
#define SIDE_COUNT_OFFSET   3

NodeArray* NodeArray_new(uint16_t capacity) {
    NodeArray *result;

    result = malloc(sizeof(*result) + capacity * sizeof(result->items[0]));
    result->length = 0;
    result->capacity = capacity;

    return result;
}

void NodeArray_free(NodeArray *array)
{
    if (array) {
        array->length = array->capacity = 0;
        free(array);
    }
}

_Bool NodeArray_is_full(NodeArray *array)
{
    return array->length == array->capacity;
}

_Bool NodeArray_is_cyclic(NodeArray *array)
{
    return 
        array->items[0]->incoming == array->items[array->length - 1]->outgoing;
}

void NodeArray_add(NodeArray *array, Node *node)
{
    assert(! NodeArray_is_full(array));
    array->items[array->length++] = node;
}

uint16_t NodeArray_sum_and_print(NodeArray *array)
{
    uint16_t result, number;
    uint8_t i;

    result = 0;
    for (i = 0; i < array->length; ++i) {
        number = array->items[i]->number;
        printf("%d ", number);
        result += number;
    }
    putchar('\n');
    return result;
}

Node* Node_new(Node *previous, uint8_t side_count, uint16_t number)
{
    Node *result;
    div_t div_result;

    result = malloc(sizeof(*result));
    result->next = previous;
    result->side_count = side_count;
    result->number = number;
    div_result = div(number, 100);
    result->incoming = div_result.quot;
    result->outgoing = div_result.rem;
    result->neighbours = NULL;

    return result;
}

void Node_free(Node *node)
{
    Node *next;

    while (node) {
        node->neighbours = NULL;
        next = node->next;
        free(node);
        node = next;
    }
}

Node* create_nodes_from_function(
    Node *node, uint8_t side_count, Function function
) {
    uint16_t n, number;

    for (n = 0;; ++n) {
        number = function(n);
        if (number >= 1000) {
            if (number >= 10000) return node;
            node = Node_new(node, side_count, number);
        }
    }
}

Node* create_nodes(const Function *functions, uint8_t side_count_offset)
{
    Node *nodes;
    uint8_t i;

    nodes = NULL;
    for (i = 0; i < FUNCTION_COUNT; ++i) {
        nodes = create_nodes_from_function(
            nodes, i + side_count_offset, functions[i]
        );
    }

    return nodes;
}

void build_graph(Node *nodes, NodeArray **incoming2nodes)
{
    Node *node;
    uint16_t incoming_histogram[100];
    uint8_t i;

    memset(incoming_histogram, 0, sizeof(incoming_histogram));
    node = nodes;
    while (node) {
        ++incoming_histogram[node->incoming];
        node = node->next;
    }

    for (i = 0; i < 100; ++i) {
        incoming2nodes[i] = NodeArray_new(incoming_histogram[i]);
    }

    node = nodes;
    while (node) {
        NodeArray_add(incoming2nodes[node->incoming], node);
        node->neighbours = incoming2nodes[node->outgoing];
        node = node->next;
    }   
}

_Bool _find_solution(Node *node, NodeArray *result, _Bool *side_counts)
{
    Node *neighbour;
    uint16_t i;

    NodeArray_add(result, node);
    side_counts[node->side_count] = true;
    if (NodeArray_is_full(result)) {
        if (NodeArray_is_cyclic(result)) {
            return true;
        }
    } else {
        for (i = 0; i < node->neighbours->length; ++i) {
            neighbour = node->neighbours->items[i];
            if (! side_counts[neighbour->side_count]) {
                if (_find_solution(neighbour, result, side_counts)) {
                    return true;
                }
            }
        }
    }
    --result->length;  /* Discard last element of `result`. */
    side_counts[node->side_count] = false;
    return false;
}

NodeArray* find_solution(Node *nodes)
{
    NodeArray *result;
    _Bool side_counts[FUNCTION_COUNT + SIDE_COUNT_OFFSET];
    Node *node;

    result = NodeArray_new(FUNCTION_COUNT);
    memset(side_counts, 0, sizeof(side_counts));
    
    node = nodes;
    while (node) {
        if (_find_solution(node, result, side_counts)) return result;
        node = node->next;
    }

    NodeArray_free(result);
    return NULL;
}

int main(void)
{
    Node *nodes;
    NodeArray *incoming2nodes[100];
    NodeArray *result;
    uint8_t i;
    
    nodes = create_nodes(FUNCTIONS, SIDE_COUNT_OFFSET);
    build_graph(nodes, incoming2nodes);
    result = find_solution(nodes);
    if (result) {
        printf("sum = %d\n", NodeArray_sum_and_print(result));
    } else {
        puts("No solution found!");
    }

    NodeArray_free(result);
    for (i = 0; i < 100; ++i) {
        NodeArray_free(incoming2nodes[i]);
    }
    Node_free(nodes);
    return 0;
}