Skip to content

Instantly share code, notes, and snippets.

@leptos-null
Last active September 21, 2018 15:23
Show Gist options
  • Save leptos-null/ceda36ceafa6d19b75aedf1775954ca1 to your computer and use it in GitHub Desktop.
Save leptos-null/ceda36ceafa6d19b75aedf1775954ca1 to your computer and use it in GitHub Desktop.
Implementation of Descartes' Rule of Signs and Rational Zeros Theorem to show information about a given rational polynomial
//
// 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));
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment