Implementation of Descartes' Rule of Signs and Rational Zeros Theorem to show information about a given rational polynomial

polymaths.c

//
//  main.c
//  polymaths
//
//  Created by Leptos on 9/18/18.
//  Copyright © 2018 Leptos. All rights reserved.
//

#include <stdio.h>
#include <math.h>
#include <float.h>
#include <stdbool.h>
#include <assert.h>

#ifndef UsePolynomialEquationPrinting
#define UsePolynomialEquationPrinting 1
#endif

/// MARK: Polynomial description
/// An array of signed integers representing the coefficients of a polynomial
/// The coefficients must be in descending order by powers
/// i.e. 3x⁴ - 5x³ + 25x² - 45x - 18 is { 3, -5, 25, -45, -18 }
/// i.e. x⁴ - 5x² + 12x is { 1, 0, -5, 12, 0 }

/**
 @brief Calculate the number of positive zeros for a given polynomial using Descartes' Rule of Signs
 @param c A polynomial as described by Polynomial description
 @param l c array length
 @return Number of possible positive zeros
 */
static uint16_t maxNumberOfPositive(int32_t *c, size_t l) {
    bool p = (c[0] > 0);
    uint16_t r = 0;
    for (typeof(l) i = 1; i < l; i++) {
        if (c[i]) {
            bool s = (c[i] > 0);
            if (p != s) {
                p ^= true;
                r++;
            }
        }
    }
    return r;
}

/**
 @brief Calculate the number of negative zeros for a given polynomial using Descartes' Rule of Signs
 @param c A polynomial as described by Polynomial description
 @param l c array length
 @return Number of possible negative zeros
 */
static uint8_t maxNumberOfNegative(int32_t *c, size_t l) {
    int32_t cc[l];
    const bool startIsEven = l%2;
    for (typeof(l) i = 0; i < l; i++) {
        cc[i] = (i%2 && startIsEven) ? c[i] : -c[i];
    }
    
    return maxNumberOfPositive(cc, l);
}

#if UsePolynomialEquationPrinting
/**
 @brief Print the standard representation of a given polynomial
 @param c A polynomial as described by Polynomial description
 @param l c array length
 */
static void putEquation(int32_t *c, size_t l) {
    static const char *superScripts[10] = { "⁰", "¹", "²", "³", "⁴", "⁵", "⁶", "⁷", "⁸", "⁹" };
    assert(l < 10 &&  "putEquation currently does not support exponents more than 1 decimal digit");
    typeof(l) u = l;
    for (typeof(l) i = 0; i < l; i++) {
        u--;
        if (c[i]) {
            bool isPositive = c[i] > 0;
            printf("%s%d%s%s", (u == l-1) ? "" : (isPositive ? " + " : " - "),
                   isPositive ? c[i] : -c[i], u ? "x" : "", (!u || u == 1) ? "" : superScripts[u]);
        }
    }
    putchar('\n');
}
#endif

/**
 @brief Get a list of factors for a given integer
 @param n The number to get factors
 @param l Buffer to copy factor list into
 @return The amount of factors
 */
static size_t getFactors(uint32_t n, uint32_t *l) {
    size_t r = 0;
    for (typeof(n) i = 1; i <= n; i++) {
        if (n%i == 0) {
            if (l) {
                l[r] = i;
            }
            r++;
        }
    }
    return r;
}

/**
 @brief Calculate the result of a polynomial given an x value
 @param c A polynomial as described by Polynomial description
 @param l c array length
 @param x x value
 @return computed result
 */
static double checkAttempt(int32_t *c, size_t l, double x) {
    double r = 0;
    typeof(l) u = l;
    for (typeof(l) i = 0; i < l; i++) {
        r += (c[--u] * pow(x, i));
    }
    return r;
}

/**
 @brief Print the integer or fractional representation of the n/d if it's a zero of a given polynomial
 @param c A polynomial as described by Polynomial description
 @param l c array length
 @param n numerator
 @param d denominator
 @param s pointer to a value that will be set to whether the value resulted in a zero
 @return x value used to calculate
 */
static double printCheckWithValues(int32_t *c, size_t l, double n, uint32_t d, bool *s) {
    if (fabs(checkAttempt(c, l, n/d)) < DBL_EPSILON) {
        if (s) {
            *s = true;
        }
        if (d == 1) {
            printf("zero: %+.0f\n", n);
        } else {
            printf("zero: %+.0f/%u\n", n, d);
        }
    } else {
        if (s) {
            *s = false;
        }
    }
    return n/d;
}

/**
 @brief Print information about a polynomial function
 @param c A polynomial as described by Polynomial description
 @param l c array length
 */
void printPlyInfo(int32_t *c, size_t l) {
    assert(c[0] && "Leading coefficient must be non-zero");
#if UsePolynomialEquationPrinting
    printf("Equation: ");
    putEquation(c, l);
#endif
    
    printf("Max positive zeros: %u\n", maxNumberOfPositive(c, l));
    printf("Max negative zeros: %u\n", maxNumberOfNegative(c, l));
    
    int32_t nv = c[l-1];
    if (nv < 0) {
        nv = -nv;
    }
    size_t nfl = getFactors(nv, NULL);
    uint32_t nf[nfl];
    getFactors(nv, nf);
    
    int32_t dv = c[0];
    if (dv < 0) {
        dv = -dv;
    }
    size_t dfl = getFactors(dv, NULL);
    uint32_t df[dfl];
    getFactors(dv, df);
    /* Rational Zeros Theorem */
    for (typeof(nfl) nfi = 0; nfi < nfl; nfi++) {
        for (typeof(dfl) dfi = 0; dfi < dfl; dfi++) {
            double nfv = nf[nfi];
            uint32_t dfv = df[dfi];
            /* TODO: check success and maintain array to avoid duplicate answers */
            printCheckWithValues(c, l,  nfv, dfv, NULL);
            printCheckWithValues(c, l, -nfv, dfv, NULL);
        }
    }
}

int main() {
    // int32_t c[] = { 1, 6, -11, -24, 28 };
    // int32_t c[] = { 1, -3, 15, -75, -250 };
    // int32_t c[] = { 2, 1, -33, -26, 56 };
    int32_t c[] = { 3, -5, 25, -45, -18 };
    
    printPlyInfo(c, sizeof(c)/sizeof(int32_t));
}