New polynomial module (not quite CCAN format).

2 files changed, 562 insertions, 0 deletions

diff --git a/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt.h b/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt.h
new file mode 100644
index 00000000..97c42830
--- /dev/null
+++ b/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt.h
@@ -0,0 +1,430 @@
+/* Polynomial ADT 

+** A polynomial module with

+** ability to add,sub,mul derivate/integrate, compose ... polynomials 

+** ..expansion in progress ...

+ * Copyright (c) 2009 I. Soule

+ * All rights reserved.

+ *

+ * Redistribution and use in source and binary forms, with or without

+ * modification, are permitted provided that the following conditions

+ * are met:

+ * 1. Redistributions of source code must retain the above copyright

+ *    notice, this list of conditions and the following disclaimer.

+ * 2. Redistributions in binary form must reproduce the above copyright

+ *    notice, this list of conditions and the following disclaimer in the

+ *    documentation and/or other materials provided with the distribution.

+ *

+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND

+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE

+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE

+ * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE

+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL

+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS

+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)

+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT

+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY

+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF

+ * SUCH DAMAGE.

+ *

+**          iasoule32@gmail.com

+*/

+

+#ifndef __POLYNOMIAL_ADT

+#define __POLYNOMIAL_ADT

+

+#include <assert.h>

+#include <stdlib.h>

+#include <stdbool.h> //C99 compliance

+#include <math.h>

+

+#define max(a, b) (a) > (b) ? (a) : (b)

+#define sgn(a)    (a) < 0 ? '+' : '-' //for quadratic factored form

+

+typedef struct node {

+    int exp;

+    float coeff;

+    struct node *next;

+}Node;

+

+typedef struct polynomial_adt {

+    Node *head;

+    int terms, hp; //hp highest power

+}PolyAdt;

+

+void display_poly(const PolyAdt *a);

+/**

+* create_adt - create a polynomial on the heap

+* @hp: the highest power in the polynomial

+*/

+PolyAdt *create_adt(int hp) 

+{

+    PolyAdt *pAdt = malloc(sizeof(PolyAdt));

+    assert(pAdt != NULL);

+    

+    pAdt->head = NULL;

+    pAdt->terms = 0;

+    pAdt->hp = hp;

+

+    return pAdt;

+}

+/**

+* create_node - creates a Node (exponent, constant and next pointer) on the heap 

+* @constant: the contant in the term 

+* @exp:      the exponent on the term

+* @next:     the next pointer to another term in the polynomial

+* 

+* This should not be called by client code (hence marked static)

+* used to assist insert_term()

+*/ 

+static Node *create_node(float constant, int exp, Node *next)

+{

+    Node *nNode = malloc(sizeof(Node));

+    assert(nNode != NULL);

+    

+    nNode->exp = exp;

+    nNode->coeff = constant;

+    nNode->next = next;

+    return nNode;

+}

+

+/**

+* insert_term - inserts a term into the polynomial

+* @pAdt: the polynomial 

+* @c:    constant value on the term

+* @e:    the exponent on the term 

+*/

+void insert_term(PolyAdt *pAdt, float c, int e)

+{

+    assert(pAdt != NULL); //assume client code didnt call create_adt()

+    Node *n = malloc(sizeof(Node));

+       

+    if(pAdt->head == NULL)

+        pAdt->head = create_node(c, e, pAdt->head);

+    else

+        for(n = pAdt->head; n->next != NULL; n = n->next); //go to the end of list

+            n->next = create_node(c, e, NULL);

+    

+    pAdt->terms++;

+}

+

+/**

+* polyImage - returns an image (direct) copy of the polynomial

+* @orig: the polynomial to be duplicated

+*/

+PolyAdt *polyImage(const PolyAdt *orig)

+{

+    PolyAdt *img = create_adt(orig->hp);

+    Node *origHead = orig->head;

+    

+    for(; origHead; origHead = origHead->next)

+             insert_term(img, origHead->coeff, origHead->exp);

+    return img;

+}

+

+

+/**

+* add - adds two polynomials together, and returns their sum (as a polynomial)

+* @a: the 1st polynomial 

+* @b: the 2nd polynomial

+*/

+PolyAdt *add(const PolyAdt *a, const PolyAdt *b)

+{

+    PolyAdt *sum;

+    Node *n, *np;

+    _Bool state = true;

+    

+    assert(a != NULL && b != NULL);

+    

+    int hpow = max(a->hp, b->hp);

+    sum = create_adt(hpow); //create space for it

+    

+    /* using state machine to compare the poly with the most terms to 

+    ** the poly with fewer, round robin type of effect comparison of

+    ** exponents => 3 Cases: Equal, Less, Greater

+    */

+        n = a->head; np = b->head;

+        while(state) {

+            /* compare the exponents */

+            if(n->exp == np->exp){

+                insert_term(sum, n->coeff + np->coeff, n->exp);

+                n = n->next;

+                np = np->next;

+            }

+            

+            else if(n->exp < np->exp){

+                insert_term(sum, np->coeff, np->exp);

+                np = np->next; //move to next term of b

+            }

+            

+            else { //greater than

+                insert_term(sum, n->coeff, n->exp);

+                n = n->next;

+            }

+            /* check whether at the end of one list or the other */

+            if(np == NULL && state == true){ //copy rest of a to sum

+                for(; n != NULL; n = n->next)

+                    insert_term(sum, n->coeff, n->exp);

+                state = false;

+            }

+            

+           if(n == NULL && state == true){

+                for(; np != NULL; np = np->next)

+                    insert_term(sum, np->coeff, np->exp);

+                state = false;

+            }       

+     }        

+    return sum;               

+}            

+

+/**

+* sub - subtracts two polynomials, and returns their difference (as a polynomial)

+* @a: the 1st polynomial 

+* @b: the 2nd polynomial

+* Aids in code reuse by negating the terms (b) and then calls the add() function

+*/

+PolyAdt *subtract(const PolyAdt *a, const PolyAdt *b)

+{

+	assert(a != NULL && b != NULL);

+

+    PolyAdt *tmp = create_adt(b->hp);

+    Node *bptr;

+    

+    for(bptr = b->head; bptr != NULL; bptr = bptr->next)

+        insert_term(tmp,-bptr->coeff,bptr->exp);  //negating b's coeffs

+    return add(a,tmp);

+}

+

+/**

+* multiply - multiply two polynomials, and returns their product (as a polynomial)

+* @a: the 1st polynomial 

+* @b: the 2nd polynomial

+*/

+PolyAdt *multiply(const PolyAdt *a, const PolyAdt *b)

+{

+	assert(a != NULL && b != NULL);

+

+    //the polys are inserted in order for now

+    PolyAdt *prod = create_adt(a->head->exp + b->head->exp);

+    Node *n = a->head, *np = b->head;

+    Node *t = b->head; 

+    

+    if(a->terms < b->terms){

+        n = b->head;

+        np = t = a->head;

+    }

+    

+    for(; n != NULL; n = n->next){

+        np = t; //reset to the beginning

+        for(; np != NULL; np = np->next){ //always the least term in this loop

+                insert_term(prod, n->coeff * np->coeff, n->exp + np->exp);

+        }

+    }

+

+    return prod;       

+}

+

+/**

+* derivative - computes the derivative of a polynomial and returns the result

+* @a: the polynomial to take the derivative upon

+*/

+PolyAdt *derivative(const PolyAdt *a)

+{

+	assert(a != NULL);

+	

+	PolyAdt *deriv = create_adt(a->head->exp - 1);

+	Node *n = a->head;

+

+	for(; n != NULL; n = n->next){

+		if(n->exp == 0) break;

+		insert_term(deriv, n->coeff * n->exp, n->exp-1);

+	}

+	return deriv;

+}

+/**

+* integrate - computes the integral of a polynomial and returns the result

+* @a: the polynomial to take the integral of

+* 

+* Will compute an indefinite integral over a

+*/

+PolyAdt *integrate(const PolyAdt *a)

+{

+	assert(a != NULL);

+	

+	PolyAdt *integrand = create_adt(a->head->exp + 1);

+	Node *n;

+

+	for(n = a->head; n != NULL; n = n->next) //very simple term by term

+        insert_term(integrand, (float)n->coeff/(n->exp+1.0F), n->exp + 1);

+    

+	return integrand;

+}

+/**

+* quadratic_roots - finds the roots of the polynomial ax^2+bx+c, a != 0 && b != 0

+* @a: the polynomial

+* @real: a pointer to float of the real(R) part of a

+* @cplx: a pointer to float of the imaginary(I) part of a

+*

+* Usage:

+* Two options can be done by the client 

+* 1. Either pass NULL to real and cplx

+*    this will display the roots by printf

+*    quadratic_roots(myPolynomial, NULL, NULL);

+*

+* 2. Pass in pointers** to type float of the real and complex

+*    if the discriminant is >0 cplx = -ve root of X

+*    quadratic_roots(myPolynomial, &realPart, &complexPart);

+*/

+void quadratic_roots(const PolyAdt *a, float *real, float *cplx)

+{

+	assert(a != NULL);

+	

+	float dscrmnt, _a, b, c;

+	float u, v;

+    

+    Node *n = a->head;

+    _a = n->coeff; b = n->next->coeff; c = n->next->next->coeff;

+    

+	dscrmnt = (b*b) - 4*_a*c;

+    u = -b/(2*_a); v = sqrt((double)fabs(dscrmnt))/(2*_a);

+    

+	if(real && !cplx || !real && cplx)

+		assert(true);

+

+	if(real == NULL && cplx == NULL){

+		if(a->hp != 2 && a->terms < 3){

+		  printf("Invalid Quadratic*, A and B must be non-zero");

+			return;

+        }

+        

+		if(dscrmnt != 0)

+			printf("X = %.2f +/- %.2f%c\n",u,v,dscrmnt < 0 ? 'I':' ');

+		else{

+			printf("(X %c %.2f)(X %c %.2f)\n",sgn(u),fabs(u),sgn(u),fabs(u));

+			printf("X1,2 = %.2f\n",u);

+		}

+	}

+	//save values in pointers

+	else {

+		if(dscrmnt < 0){ //x = u +/- vI Re(x) = u, Im(x) = +v

+			*real = u; 

+			*cplx = v; //understand +/- is not representable

+		}

+		else if(dscrmnt == 0){

+			*real = u; 

+			*cplx = 0.00F;

+		}

+		else{

+			*real = u + v;

+			*cplx = u - v;

+		}

+	}

+}

+

+/**

+* exponentiate - computes polynomial exponentiation (P(x))^n, n E Z*

+* @a: the polynomial

+* @n: the exponent

+* Works fast for small n (n < 8) currently runs ~ O(n^2 lg n)

+*/

+PolyAdt *exponentiate(const PolyAdt *a, int n)

+{

+	assert(a != NULL);

+

+	PolyAdt *expn = create_adt(a->hp *  n);

+	PolyAdt *aptr = polyImage(a);

+    int hl = n / 2;

+    

+    //check default cases before calculation

+    if(n == 0){

+        insert_term(expn, 1, 0);

+        return expn;

+    }

+    else if(n == 1){

+        return aptr;

+    }

+        

+	for(; hl ; hl--)

+        aptr = multiply(aptr, aptr);

+

+    if(n % 2) //odd exponent do a^(n-1) * a = a^n

+        expn = multiply(aptr, a);

+    else

+        expn = aptr;

+    return expn;

+}

+/**

+* compose - computes the composition of two polynomials P(Q(x)) and returns the composition

+* @p: polynomial P(x) which will x will be equal to Q(x)

+* @q: polynomial Q(x) which is the argument to P(x)

+*/

+PolyAdt *compose(const PolyAdt *p, const PolyAdt *q)

+{

+    assert(p && q);

+    

+	PolyAdt *comp = create_adt(p->head->exp * q->head->exp);

+	PolyAdt *exp;

+	

+    Node *pp = p->head;

+    Node *qq = q->head;

+    

+    int swap = 0;

+    

+    if(p->terms < q->terms){

+        pp = q->head;

+        qq = p->head;

+        swap = 1;

+    }

+    

+    /* going through, exponentiate each term with the exponent of p */

+        for(; pp != NULL; pp = pp->next){

+                exp = exponentiate(swap ? p: q, pp->exp);

+                insert_term(comp, pp->coeff * exp->head->coeff, exp->head->exp);

+        }

+    

+    return comp;

+}

+/** 

+* destroy_poly - completely frees the polynomial from the heap and resets all values

+* @poly: the polynomial to release memory back to the heap

+* Usage:

+* destroy_poly(myPoly); //puts polynomial on free list

+*/

+void destroy_poly(PolyAdt *poly)

+{

+    Node *ps = poly->head;

+    Node *tmp = NULL;

+    while(ps != NULL){

+        tmp = ps;

+        free(tmp);

+        ps = ps->next;

+    }

+    poly->hp = poly->terms = 0;

+    poly->head = NULL;

+}

+/**

+* display_poly - displays the polynomial to the console in nice format

+* @a: the polynomial to display 

+*/

+void display_poly(const PolyAdt *a)

+{

+    assert(a != NULL);

+    Node *n;

+    

+    for(n = a->head; n != NULL; n = n->next){

+        

+       n->coeff < 0 ? putchar('-') : putchar('+'); 

+        if(n->exp == 0)

+            printf(" %.2f ",fabs(n->coeff));

+        else if(n->coeff == 1)

+            printf(" X^%d ",n->exp);

+        else if(n->exp == 1)

+            printf(" %.2fX ",fabs(n->coeff));

+        else if(n->coeff == 0)

+            continue;

+        else

+            printf(" %.2fX^%d ",fabs(n->coeff),n->exp);

+        }

+    printf("\n\n");

+}

+

+#endif

diff --git a/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt_test.c b/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt_test.c
new file mode 100644
index 00000000..f7c9c8a5
--- /dev/null
+++ b/junkcode/iasoule32@gmail.com-polynomial/polynomial_adt_test.c
@@ -0,0 +1,132 @@
+/* 

+** polynomial_adt_test.c

+** Test (minimalistic) for the polynomial module

+ * More of a display of functionality

+ * Copyright (c) 2009 I. Soule

+ * All rights reserved.

+ *

+ * Redistribution and use in source and binary forms, with or without

+ * modification, are permitted provided that the following conditions

+ * are met:

+ * 1. Redistributions of source code must retain the above copyright

+ *    notice, this list of conditions and the following disclaimer.

+ * 2. Redistributions in binary form must reproduce the above copyright

+ *    notice, this list of conditions and the following disclaimer in the

+ *    documentation and/or other materials provided with the distribution.

+ *

+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND

+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE

+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE

+ * ARE DISCLAIMED.  IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE

+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL

+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS

+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)

+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT

+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY

+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF

+ * SUCH DAMAGE.

+ *

+**          iasoule32@gmail.com

+*/

+

+#include <stdio.h>

+#include "polynomial_adt.h"

+

+int main()

+{

+    PolyAdt *p = create_adt(5), *q = create_adt(4);

+    PolyAdt *sum, *diff, *prod;

+    

+    insert_term(p,1,5);

+    insert_term(p,3,4);

+    insert_term(p,1,3);

+    insert_term(p,9,2);

+    insert_term(p,8,1);

+    

+    insert_term(q,2,4);

+    insert_term(q,8,3);

+    insert_term(q,7,2);

+    insert_term(q,6,1);

+    

+    

+    printf("Displaying Polynomials ...\n");

+    display_poly(p);

+    display_poly(q);

+    

+    sum = add(p,q);

+    printf("P(x) + Q(x) = ");

+    display_poly(sum);

+    

+    diff = subtract(p,q);

+    printf("P(x) - Q(x) = ");

+    display_poly(diff);

+    

+    prod = multiply(p,q);

+    printf("P(x)*Q(x) = ");

+    display_poly(prod);

+    

+    PolyAdt *quad = create_adt(2);

+    insert_term(quad, 10, 2);

+    insert_term(quad, 30, 1);

+    insert_term(quad, 2, 0);

+    

+    quadratic_roots(quad, NULL, NULL); //print out the roots

+    

+    float real, cplx;

+    quadratic_roots(quad, &real, &cplx);

+    

+    printf("X1 = %f, X2 = %f\n\n", real, cplx);

+    

+    PolyAdt *deriv, *integral;

+    

+    deriv = derivative(p);

+    printf("The derivitive of p = ");

+    display_poly(deriv);

+    integral = integrate(q);

+    

+    printf("The integral of q = ");

+    display_poly(integral);

+    

+    printf("\n Computing P(x)^3\n");

+    

+    PolyAdt *expo;

+    expo = exponentiate(p, 3);

+    display_poly(expo);

+    printf("\n");

+    

+    printf("Computing Integral[Q(x)^2]\n");

+    expo = exponentiate(q, 2);

+    integral = integrate(expo);

+    display_poly(integral);

+    

+    

+    printf(" Differentiating and Integrating P\n");

+    display_poly(integrate(derivative(p))); 

+    

+    PolyAdt *L, *M;

+    

+    L = create_adt(3), M = create_adt(2);

+    

+    insert_term(L, 4, 3);

+    insert_term(L, 10, 2);

+    insert_term(L, 15, 1);

+    

+    insert_term(M, 4, 2);

+    printf("L = ");

+    display_poly(L);

+    printf("M = ");

+    display_poly(M);

+    

+    

+    printf("Computing composition L(M(X))\n");

+    display_poly(compose(L, M));

+    

+    printf("Freed memory back to heap for allocated polys'");

+    destroy_poly(sum);

+    destroy_poly(diff);

+    destroy_poly(prod);

+    destroy_poly(L); destroy_poly(M);

+    destroy_poly(q); destroy_poly(p);

+    

+    return 0;

+}