Initial commit.

1 files changed, 551 insertions, 0 deletions

diff --git a/ic-reals-6.3/base/gt0.c b/ic-reals-6.3/base/gt0.c
new file mode 100644
index 0000000..83ed0af
--- /dev/null
+++ b/ic-reals-6.3/base/gt0.c
@@ -0,0 +1,551 @@
+/*
+ * Copyright (C) 2000, Imperial College
+ *
+ * This file is part of the Imperial College Exact Real Arithmetic Library.
+ * See the copyright notice included in the distribution for conditions
+ * of use.
+ */
+
+#include <stdio.h>
+#include "real.h"
+#include "real-impl.h"
+
+/*
+ * This file contains everything to do with the predicate x > 0.
+ * This can be defined in terms of x >= 0 since
+ * x > 0 \equiv not_B(-x >= 0) but no matter.
+ */
+
+void setPredX(PredX *, BoolVal);
+void absorbDigsXIntoPredX(PredX *);
+void absorbSignXIntoPredX(PredX *);
+void force_To_PredX_From_The_Abyss();
+static BoolVal gt_Vec_0(Vec *);
+static BoolVal gt_MatX_0(MatX *);
+
+#define forceGtZero_To_PredX_From_DigsX_SINF_Entry \
+			forceGtEqZero_To_PredX_From_DigsX_SINF_Entry
+#define forceGtZero_To_PredX_From_DigsX_SNEG_Entry \
+			force_To_PredX_From_DigsX_minus_2n_minus_1_False_Entry
+
+void setGtZeroPredXMethod(PredX *);
+
+/*
+ * You'll find slightly more comments in the file gteq0.c
+ */
+Bool
+gt_R_0(Real x)
+{
+	PredX *predX;
+	PredX *allocPredX();
+	void setGtZeroPredXMethod(PredX *);
+
+	predX = allocPredX(x);
+	setGtZeroPredXMethod(predX);
+	return (Bool) predX;
+}
+
+void
+setGtZeroPredXMethod(PredX *predX)
+{
+	void forceGtZero_To_PredX_From_SignX_Entry();
+	void forceGtZero_To_PredX_From_DigsX_SPOS_Entry();
+	void forceGtZero_To_PredX_From_MatX_Entry();
+	void forceGtZero_To_PredX_From_MatX_Signed_Entry();
+	void forceGtZero_To_PredX_From_TenXY();
+	void forceGtZero_To_PredX_From_Alt_Entry();
+	void forceGtZero_To_PredX_From_Cls_Entry();
+
+	switch (predX->x->gen.tag.type) {
+	case VECTOR :
+		setPredX(predX, gt_Vec_0((Vec *) predX->x));
+		predX->force = force_To_PredX_From_The_Abyss;
+		break;
+	case MATX :
+		if (predX->x->matX.x->gen.tag.isSigned)
+			predX->force = forceGtZero_To_PredX_From_MatX_Signed_Entry;
+		else {
+			setPredX(predX, gt_MatX_0((MatX *) predX->x));
+			predX->force = forceGtZero_To_PredX_From_MatX_Entry;
+		}
+		break;
+	case TENXY :
+		predX->force = forceGtZero_To_PredX_From_TenXY;
+		break;
+	case SIGNX :
+		predX->force = forceGtZero_To_PredX_From_SignX_Entry;
+		break;
+	case DIGSX :
+		predX->force = forceGtZero_To_PredX_From_DigsX_SPOS_Entry;
+		break;
+	case CLOSURE :
+		predX->force = forceGtZero_To_PredX_From_Cls_Entry;
+		break;
+	case ALT :
+		predX->force = forceGtZero_To_PredX_From_Alt_Entry;
+		break;
+	default :
+		Error(FATAL, E_INT, "compareGtZero", "argument is not a stream");
+	}
+}
+
+void
+forceGtZero_To_PredX_From_Alt_Entry()
+{
+	PredX *predX;
+	Alt *alt;
+	void force_To_Alt_Entry();
+	void forceGtZero_To_PredX_From_Alt_Cont();
+
+	predX = (PredX *) POP;
+	alt = (Alt *) predX->x;
+	
+	PUSH_2(forceGtZero_To_PredX_From_Alt_Cont, predX);
+
+	/*
+	 * If alt->redirect is not valid (equals NULL) then the value of
+	 * the conditional has not been determined so we need to force it.
+	 */
+	if (alt->redirect == NULL)
+		PUSH_2(force_To_Alt_Entry, alt);
+}
+
+void
+forceGtZero_To_PredX_From_Alt_Cont()
+{
+	PredX *predX;
+	Alt *alt;
+
+	predX = (PredX *) POP;
+	alt = (Alt *) predX->x;
+
+#ifdef DAVINCI
+	beginGraphUpdate();
+	deleteOnlyEdge(predX, alt);
+	newEdgeToOnlyChild(predX, alt->redirect);
+	endGraphUpdate();
+#endif
+	predX->x = alt->redirect;
+	setGtZeroPredXMethod(predX);
+
+	PUSH_2(predX->force, predX);
+}
+
+void
+forceGtZero_To_PredX_From_Cls_Entry()
+{
+	PredX *predX;
+	Cls *cls;
+	void forceGtZero_To_PredX_From_Cls_Cont();
+
+	predX = (PredX *) POP;
+	cls = (Cls *) predX->x;
+	
+	PUSH_2(forceGtZero_To_PredX_From_Cls_Cont, predX);
+
+	/*
+	 * If cls->redirect is not valid (equals NULL) then the value of
+	 * the closure has not been determined so we need to force it.
+	 */
+	if (cls->redirect == NULL)
+		PUSH_2(cls->force, cls);
+}
+
+void
+forceGtZero_To_PredX_From_Cls_Cont()
+{
+	PredX *predX;
+	Cls *cls;
+
+	predX = (PredX *) POP;
+	cls = (Cls *) predX->x;
+
+#ifdef DAVINCI
+	beginGraphUpdate();
+	deleteOnlyEdge(predX, cls);
+	newEdgeToOnlyChild(predX, cls->redirect);
+	endGraphUpdate();
+#endif
+	predX->x = cls->redirect;
+	setGtZeroPredXMethod(predX);
+
+	PUSH_2(predX->force, predX);
+}
+
+/*
+ * If we get here then we need to get more information into the matrix.
+ */
+void
+forceGtZero_To_PredX_From_MatX_Entry()
+{
+	PredX *predX;
+	MatX *matX;
+	void forceGtZero_To_PredX_From_MatX_Cont();
+
+	predX = (PredX *) POP;
+	matX = (MatX *) predX->x;
+
+	if (matX->tag.type == VECTOR) {
+		setPredX(predX, gt_Vec_0((Vec *) matX));
+		predX->force = force_To_PredX_From_The_Abyss;
+		return;
+	}
+
+	PUSH_2(forceGtZero_To_PredX_From_MatX_Cont, predX);
+	PUSH_3(predX->x->matX.force, predX->x, defaultForceCount);
+}
+
+void
+forceGtZero_To_PredX_From_MatX_Cont()
+{
+	PredX *predX;
+	MatX *matX;
+
+	predX = (PredX *) POP;
+	matX = (MatX *) predX->x;
+
+	if (matX->tag.type == VECTOR) {
+		setPredX(predX, gt_Vec_0((Vec *) matX));
+		predX->force = force_To_PredX_From_The_Abyss;
+		return;
+	}
+
+	setPredX(predX, gt_MatX_0((MatX *) (predX->x)));
+}
+
+void
+forceGtZero_To_PredX_From_MatX_Signed_Entry()
+{
+	PredX *predX;
+	MatX *matX;
+	void forceGtZero_To_PredX_From_MatX_Signed_Cont();
+
+	predX = (PredX *) POP;
+	matX = (MatX *) predX->x;
+
+	if (matX->tag.type == VECTOR) {
+		setPredX(predX, gt_Vec_0((Vec *) matX));
+		predX->force = force_To_PredX_From_The_Abyss;
+		return;
+	}
+
+	PUSH_2(forceGtZero_To_PredX_From_MatX_Signed_Cont, predX);
+
+	if (matX->x->gen.tag.isSigned)
+		PUSH_2(matX->force, matX);
+}
+
+void
+forceGtZero_To_PredX_From_MatX_Signed_Cont()
+{
+	PredX *predX;
+	MatX *matX;
+	void forceGtZero_To_PredX_From_MatX_Entry();
+
+	predX = (PredX *) POP;
+	matX = (MatX *) predX->x;
+
+	if (matX->tag.type == VECTOR) {
+		setPredX(predX, gt_Vec_0((Vec *) matX));
+		predX->force = force_To_PredX_From_The_Abyss;
+		return;
+	}
+	predX->force = forceGtZero_To_PredX_From_MatX_Entry;
+	setPredX(predX, gt_MatX_0(matX));
+}
+
+void
+forceGtZero_To_PredX_From_TenXY()
+{
+	PredX *predX;
+	TenXY *tenXY;
+	void forceGtZero_To_PredX_From_SignX_Entry();
+	void forceGtZero_To_PredX_From_DigsX_SPOS_Entry();
+
+	predX = (PredX *) POP;
+	tenXY = (TenXY *) predX->x;
+
+	/*
+	 * The tensor may have reduced to a vector
+	 */
+	if (tenXY->tag.type == VECTOR) {
+		setPredX(predX, gt_Vec_0((Vec *) tenXY));
+		predX->force = force_To_PredX_From_The_Abyss;
+		return;
+	}
+
+	/*
+	 * The tensor may have reduced to a matrix (signed or otherwise)
+	 */
+	if (tenXY->tag.type == MATX) {
+		if (predX->x->matX.x->gen.tag.isSigned) {
+			predX->force = forceGtZero_To_PredX_From_MatX_Signed_Entry;
+			PUSH_2(predX->force, predX);
+		}
+		else {
+			setPredX(predX, gt_MatX_0((MatX *) predX->x));
+			predX->force = forceGtZero_To_PredX_From_MatX_Entry;
+		}
+		return;
+	}
+
+	if (tenXY->tag.isSigned) {
+		createSignedStreamForTenXY(tenXY);
+#ifdef DAVINCI
+		beginGraphUpdate();
+		deleteOnlyEdge(predX, predX->x);
+		newEdgeToOnlyChild(predX, tenXY->strm);
+		endGraphUpdate();
+#endif
+		predX->x = tenXY->strm;
+		predX->force = forceGtZero_To_PredX_From_SignX_Entry;
+		PUSH_2(predX->force, predX);
+	}
+	else {
+		createUnsignedStreamForTenXY(tenXY);
+#ifdef DAVINCI
+		beginGraphUpdate();
+		deleteOnlyEdge(predX, predX->x);
+		newEdgeToOnlyChild(predX, tenXY->strm);
+		endGraphUpdate();
+#endif
+		predX->x = tenXY->strm;
+		predX->force = forceGtZero_To_PredX_From_DigsX_SPOS_Entry;
+		PUSH_2(predX->force, predX);
+	}
+}
+
+void
+forceGtZero_To_PredX_From_SignX_Entry()
+{
+	PredX *predX;
+	SignX *signX;
+	void forceGtZero_To_PredX_From_SignX_Cont();
+
+	predX = (PredX *) POP;
+	signX = (SignX *) predX->x;
+
+	PUSH_2(forceGtZero_To_PredX_From_SignX_Cont, predX);
+
+	if (signX->tag.value == SIGN_UNKN)
+		PUSH_2(signX->force, signX);
+}
+
+/*
+ * At this point we know the sign is valid so we ``absorb'' it. In this
+ * case that means that we have to inspect the sign and reset the methods
+ * before advancing to the first characteristic pair.
+ */
+void
+forceGtZero_To_PredX_From_SignX_Cont()
+{
+	PredX *predX;
+	SignX *signX;
+	void forceGtZero_To_PredX_From_DigsX_SPOS_Entry();
+	void forceGtZero_To_PredX_From_DigsX_SNEG_Entry();
+	void forceGtZero_To_PredX_From_DigsX_SZERO_Entry();
+	void forceGtZero_To_PredX_From_DigsX_SINF_Entry();
+
+	predX = (PredX *) POP;
+	signX = (SignX *) predX->x;
+
+	switch (signX->tag.value) {
+	case SPOS :
+		predX->force = forceGtZero_To_PredX_From_DigsX_SPOS_Entry;
+		break;
+	case SNEG :
+		predX->force = forceGtZero_To_PredX_From_DigsX_SNEG_Entry;
+		break;
+	case SZERO :
+		predX->force = forceGtZero_To_PredX_From_DigsX_SZERO_Entry;
+		break;
+	case SINF :
+		predX->force = forceGtZero_To_PredX_From_DigsX_SINF_Entry;
+		break;
+	default :
+		Error(FATAL, E_INT, "forceGtZero_To_PredX_From_SignX_Cont",
+				"invalid sign");
+	}
+	absorbSignXIntoPredX(predX);
+	PUSH_2(predX->force, predX);
+}
+
+void
+forceGtZero_To_PredX_From_DigsX_SPOS_Entry()
+{
+	PredX *predX;
+	DigsX *digsX;
+	void forceGtZero_To_PredX_From_DigsX_SPOS_Cont();
+
+	predX = (PredX *) POP;
+	digsX = (DigsX *) predX->x;
+	
+	PUSH_2(forceGtZero_To_PredX_From_DigsX_SPOS_Cont, predX);
+	if (digsX->count == 0)
+		PUSH_3(digsX->force, digsX, defaultForceCount);
+}
+
+void
+forceGtZero_To_PredX_From_DigsX_SPOS_Cont()
+{
+	PredX *predX;
+	DigsX *digsX;
+	int k;
+	void force_To_PredX_From_DigsX_2n_minus_1_True_Entry();
+	void force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry();
+
+	predX = (PredX *) POP;
+	digsX = (DigsX *) predX->x;
+
+#ifdef PACK_DIGITS
+	if (digsX->count <= DIGITS_PER_WORD) {
+		k = (1 << digsX->count) - 1;
+		if (digsX->word.small == k) {
+			absorbDigsXIntoPredX(predX);
+			predX->force = force_To_PredX_From_DigsX_2n_minus_1_True_Entry;
+		}
+		else {
+ 			if (digsX->word.small == -k) {
+				absorbDigsXIntoPredX(predX);
+				predX->force =
+						force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry;
+			}
+			else
+				setPredX(predX, LAZY_TRUE);
+		}
+	}
+	else {
+#endif
+		/*
+		 * This is comparing a big word with +-(2^n - 1). It would be faster
+		 * to compare each word with 0xffffffff but this may have to
+		 * wait.	####
+		 */
+		if (mpz_sgn(digsX->word.big) >= 0) {
+			if (mpz_popcount(digsX->word.big) == digsX->count) {
+				absorbDigsXIntoPredX(predX);
+				predX->force = force_To_PredX_From_DigsX_2n_minus_1_True_Entry;
+			}
+			else
+				setPredX(predX, LAZY_TRUE);
+		}
+		else {
+			/*
+			 * We negate things here but I don't think it is necessary
+			 * since GMP uses sign and magnitude representation 
+			 * for big integers. Leave it for now.
+			 */
+			mpz_neg(digsX->word.big, digsX->word.big);
+			if (mpz_popcount(digsX->word.big) == digsX->count) {
+				absorbDigsXIntoPredX(predX);
+				predX->force =
+						force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry;
+			}
+			else
+				setPredX(predX, LAZY_TRUE);
+			mpz_neg(digsX->word.big, digsX->word.big);
+		}
+#ifdef PACK_DIGITS
+	}
+#endif
+}
+
+void
+forceGtZero_To_PredX_From_DigsX_SZERO_Entry()
+{
+	PredX *predX;
+	DigsX *digsX;
+	void forceGtZero_To_PredX_From_DigsX_SZERO_Cont();
+
+	predX = (PredX *) POP;
+	digsX = (DigsX *) predX->x;
+
+	PUSH_2(forceGtZero_To_PredX_From_DigsX_SZERO_Cont, predX);
+	if (digsX->count == 0)
+		PUSH_3(digsX->force, digsX, defaultForceCount);
+}
+
+void
+forceGtZero_To_PredX_From_DigsX_SZERO_Cont()
+{
+	PredX *predX;
+	DigsX *digsX;
+	void force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry();
+
+	predX = (PredX *) POP;
+	digsX = (DigsX *) predX->x;
+
+#ifdef PACK_DIGITS
+	if (digsX->count <= DIGITS_PER_WORD) {
+		if (digsX->word.small < 0)
+			setPredX(predX, LAZY_FALSE);
+		else {
+			if (digsX->word.small == 0)
+				absorbDigsXIntoPredX(predX);
+			else
+				if (digsX->word.small == 1) {
+					absorbDigsXIntoPredX(predX);
+					predX->force =
+						force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry;
+				}
+				else
+					setPredX(predX, LAZY_TRUE);
+		}
+	}
+	else {
+#endif
+		switch (mpz_sgn(digsX->word.big)) {
+		case 1 :
+			if (mpz_cmp_si(digsX->word.big, 1) > 0)
+				setPredX(predX, LAZY_TRUE);
+			else { /* word == 1 */
+				absorbDigsXIntoPredX(predX);
+				predX->force =
+						force_To_PredX_From_DigsX_minus_2n_minus_1_True_Entry;
+			}
+			break;
+		case 0 :
+			absorbDigsXIntoPredX(predX);
+			break;
+		case -1 :
+			setPredX(predX, LAZY_FALSE);
+			break;
+		default :
+			Error(FATAL, E_INT, "forceGtZero_To_PredX_From_DigsX_SZERO_Cont",
+					"bad value returned from mpz_sgn");
+		}
+#ifdef PACK_DIGITS
+	}
+#endif
+}
+
+static BoolVal
+gt_MatX_0(MatX *matX)
+{
+	if (mpz_sgn(matX->mat[0][1]) * mpz_sgn(matX->mat[1][1]) > 0) {
+		if (mpz_sgn(matX->mat[0][0]) * mpz_sgn(matX->mat[0][1]) > 0) {
+			if (mpz_sgn(matX->mat[1][0]) * mpz_sgn(matX->mat[1][1]) > 0)
+				return LAZY_TRUE;
+		}
+		else {
+			if (mpz_sgn(matX->mat[1][0]) * mpz_sgn(matX->mat[1][1]) <= 0)
+				return LAZY_FALSE;
+		}
+	}
+	return LAZY_UNKNOWN;
+}
+ 
+static BoolVal
+gt_Vec_0(Vec *vec)
+{
+	if (mpz_sgn(vec->vec[1]) != 0) {
+		if (mpz_sgn(vec->vec[0]) * mpz_sgn(vec->vec[1]) > 0)
+			return LAZY_TRUE;
+		else
+			return LAZY_FALSE;
+	}
+	if (mpz_sgn(vec->vec[0]) == 0)
+		Error(FATAL, E_INT, "gt_Vec_0", "malformed vector");
+ 
+	return LAZY_UNKNOWN;
+}