aboutsummaryrefslogtreecommitdiff
path: root/ic-reals-6.3/base/emitSign.c
diff options
context:
space:
mode:
Diffstat (limited to 'ic-reals-6.3/base/emitSign.c')
-rw-r--r--ic-reals-6.3/base/emitSign.c618
1 files changed, 618 insertions, 0 deletions
diff --git a/ic-reals-6.3/base/emitSign.c b/ic-reals-6.3/base/emitSign.c
new file mode 100644
index 0000000..d25aa03
--- /dev/null
+++ b/ic-reals-6.3/base/emitSign.c
@@ -0,0 +1,618 @@
+/*
+ * 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 deals with emitting signs from LFTs. The algorithm is similar
+ * to digit emission. See ``Notes on sign and digit emission'' for details.
+ */
+
+/*
+ * This version of emitSign favours SZERO over all other signs. This
+ * is because most of the analytic functions prefer SZERO.
+ */
+typedef struct TE {
+ Sign sign;
+ Comparison compare;
+ void (*residual)(Vector); /* the function to compute the residual */
+ struct TE *next; /* computes next entry in the sequence */
+} TblEntry;
+
+static bool mergeTwoLists(TblEntry **, TblEntry **, Comparison *);
+static bool mergeFourLists(TblEntry **, TblEntry **, TblEntry **, TblEntry **);
+static TblEntry *setList(Vector);
+
+static mpz_t tmp_z;
+
+/*
+ * Functions prefixed with 'r' compute the residual vectors. This is done
+ * "in-place". By this we mean that once we have decided that a sign 's'
+ * can be emitted from a vector 'v', then we overwrite 'v' with residual
+ * after 's' is emitted.
+ */
+
+static void rNoOp(Vector);
+
+#define rgg rNoOp
+#define regSPOS rNoOp
+#define rgeSPOS rNoOp
+
+static void regSNEG(Vector);
+static void regSZERO(Vector);
+static void relSPOS(Vector);
+static void relSNEG(Vector);
+static void relSZERO(Vector);
+static void rgeSNEG(Vector);
+static void rgeSINF(Vector);
+static void rleSPOS(Vector);
+static void rleSNEG(Vector);
+static void rleSINF(Vector);
+static void rggeSZERO(Vector);
+static void rggeSINF(Vector);
+static void rgggSINF(Vector);
+static void rgglSZERO(Vector);
+static void rll(Vector);
+static void rlleSZERO(Vector);
+static void rlleSINF(Vector);
+static void rllgSZERO(Vector);
+static void rlllSINF(Vector);
+static void rgl(Vector);
+static void rgleSZERO(Vector);
+static void rgleSINF(Vector);
+static void rglgSZERO(Vector);
+static void rgllSINF(Vector);
+static void rlg(Vector);
+static void rlgeSZERO(Vector);
+static void rlgeSINF(Vector);
+static void rlggSINF(Vector);
+static void rlglSZERO(Vector);
+
+static TblEntry tee4 = {SNEG, EQ, rNoOp, NULL};
+static TblEntry tee3 = {SPOS, EQ, rNoOp, &tee4};
+static TblEntry tee2 = {SINF, EQ, rNoOp, &tee3};
+static TblEntry tee1 = {SZERO, EQ, rNoOp, &tee2};
+
+static TblEntry teg3 = {SNEG, GT, regSNEG, NULL};
+static TblEntry teg2 = {SPOS, GT, regSPOS, &teg3};
+static TblEntry teg1 = {SZERO, GT, regSZERO, &teg2};
+
+static TblEntry tel3 = {SNEG, LT, relSNEG, NULL};
+static TblEntry tel2 = {SPOS, LT, relSPOS, &tel3};
+static TblEntry tel1 = {SZERO, LT, relSZERO, &tel2};
+
+static TblEntry tge3 = {SNEG, LT, rgeSNEG, NULL};
+static TblEntry tge2 = {SPOS, GT, rgeSPOS, &tge3};
+static TblEntry tge1 = {SINF, GT, rgeSINF, &tge2};
+
+static TblEntry tle3 = {SNEG, GT, rleSNEG, NULL};
+static TblEntry tle2 = {SPOS, LT, rleSPOS, &tle3};
+static TblEntry tle1 = {SINF, LT, rleSINF, &tle2};
+
+static TblEntry tgge3 = {SPOS, GT, rgg, NULL};
+static TblEntry tgge2 = {SINF, GT, rggeSINF, &tgge3};
+static TblEntry tgge1 = {SZERO, GT, rggeSZERO, &tgge2};
+
+static TblEntry tggg2 = {SPOS, GT, rgg, NULL};
+static TblEntry tggg1 = {SINF, GT, rgggSINF, &tggg2};
+
+static TblEntry tggl2 = {SPOS, GT, rgg, NULL};
+static TblEntry tggl1 = {SZERO, GT, rgglSZERO, &tggl2};
+
+static TblEntry tlle3 = {SPOS, LT, rll, NULL};
+static TblEntry tlle2 = {SINF, LT, rlleSINF, &tlle3};
+static TblEntry tlle1 = {SZERO, LT, rlleSZERO, &tlle2};
+
+static TblEntry tllg2 = {SPOS, LT, rll, NULL};
+static TblEntry tllg1 = {SZERO, LT, rllgSZERO, &tllg2};
+
+static TblEntry tlll2 = {SPOS, LT, rll, NULL};
+static TblEntry tlll1 = {SINF, LT, rlllSINF, &tlll2};
+
+static TblEntry tgle3 = {SNEG, LT, rgl, NULL};
+static TblEntry tgle2 = {SINF, GT, rgleSINF, &tgle3};
+static TblEntry tgle1 = {SZERO, LT, rgleSZERO, &tgle2};
+
+static TblEntry tglg2 = {SNEG, LT, rgl, NULL};
+static TblEntry tglg1 = {SZERO, LT, rglgSZERO, &tglg2};
+
+static TblEntry tgll2 = {SNEG, LT, rgl, NULL};
+static TblEntry tgll1 = {SINF, GT, rgllSINF, &tgll2};
+
+static TblEntry tlge3 = {SNEG, GT, rlg, NULL};
+static TblEntry tlge2 = {SINF, LT, rlgeSINF, &tlge3};
+static TblEntry tlge1 = {SZERO, GT, rlgeSZERO, &tlge2};
+
+static TblEntry tlgg2 = {SNEG, GT, rlg, NULL};
+static TblEntry tlgg1 = {SINF, LT, rlggSINF, &tlgg2};
+
+static TblEntry tlgl2 = {SNEG, GT, rlg, NULL};
+static TblEntry tlgl1 = {SZERO, GT, rlglSZERO, &tlgl2};
+
+void
+initEmitSign(void)
+{
+ mpz_init(tmp_z);
+}
+
+bool
+emitSignFromVector(Vector v, Sign *s)
+{
+ TblEntry *e;
+
+ if ((e = setList(v)) != NULL) {
+ (e->residual)(v);
+ *s = e->sign;
+#ifdef TRACE
+ if (TRACE && !validVector(v)) {
+ printf("emitSignFromVector: bad vector, sign=%d\n", *s);
+ dumpVector(v);
+ Error(FATAL, E_INT, "emitSignFromVector", "bad vector");
+ }
+#endif
+ return TRUE;
+ } else
+ return FALSE;
+}
+
+bool
+emitSignFromMatrix(Matrix m, Sign *s)
+{
+ TblEntry *e0, *e1;
+ Comparison comp; /* not used */
+
+ e0 = setList(m[0]);
+ e1 = setList(m[1]);
+
+ if (mergeTwoLists(&e0, &e1, &comp)) {
+ (e0->residual)(m[0]);
+ (e1->residual)(m[1]);
+ *s = e0->sign;
+#ifdef TRACE
+ if (TRACE && !validMatrix(m)) {
+ printf("emitSignFromMatrix: bad matrix, sign=%d\n", *s);
+ dumpMatrix(m);
+ Error(FATAL, E_INT, "emitSignFromMatrix", "bad matrix");
+ }
+#endif
+ return TRUE;
+ } else
+ return FALSE;
+}
+
+bool
+emitSignFromTensor(Tensor t, Sign *s)
+{
+ TblEntry *e0, *e1, *e2, *e3;
+
+ e0 = setList(t[0]);
+ e1 = setList(t[1]);
+ e2 = setList(t[2]);
+ e3 = setList(t[3]);
+
+ if (mergeFourLists(&e0, &e1, &e2, &e3)) {
+ (e0->residual)(t[0]);
+ (e1->residual)(t[1]);
+ (e2->residual)(t[2]);
+ (e3->residual)(t[3]);
+ *s = e0->sign;
+#ifdef TRACE
+ if (TRACE && !validTensor(t)) {
+ printf("emitSignFromTensor: bad tensor, sign=%d\n", *s);
+ dumpTensor(t);
+ Error(FATAL, E_INT, "emitSignFromTensor", "bad tensor");
+ }
+#endif
+ return TRUE;
+ } else
+ return FALSE;
+}
+
+static bool
+mergeFourLists(TblEntry **e0, TblEntry **e1, TblEntry **e2, TblEntry **e3)
+{
+ bool flag;
+ Comparison compLeft, compRight;
+
+ flag = mergeTwoLists(e0, e1, &compLeft)
+ && mergeTwoLists(e2, e3, &compRight);
+
+ while (flag) {
+ if ((*e0)->sign == (*e2)->sign) {
+ if (compLeft == compRight || compLeft == EQ || compRight == EQ)
+ return TRUE;
+ else {
+ *e0 = (*e0)->next;
+ *e1 = (*e1)->next;
+ *e2 = (*e2)->next;
+ *e3 = (*e3)->next;
+ flag = mergeTwoLists(e0, e1, &compLeft)
+ && mergeTwoLists(e2, e3, &compRight);
+ }
+ }
+ else {
+ if ((*e0)->sign < (*e2)->sign) {
+ *e0 = (*e0)->next;
+ *e1 = (*e1)->next;
+ flag = mergeTwoLists(e0, e1, &compLeft);
+ }
+ else {
+ *e2 = (*e2)->next;
+ *e3 = (*e3)->next;
+ flag = mergeTwoLists(e2, e3, &compRight);
+ }
+ }
+ }
+ return FALSE;
+}
+
+static bool
+mergeTwoLists(TblEntry **e0, TblEntry **e1, Comparison *comp)
+{
+ while (*e0 != NULL && *e1 != NULL) {
+ if ((*e0)->sign == (*e1)->sign) {
+ if ((*e0)->compare == (*e1)->compare) {
+ *comp = (*e0)->compare;
+ return TRUE;
+ }
+ else {
+ if ((*e0)->compare == EQ) {
+ *comp = (*e1)->compare; /* must be GT or LT */
+ return TRUE;
+ }
+ else {
+ if ((*e1)->compare == EQ) {
+ *comp = (*e0)->compare; /* must be GT or LT */
+ return TRUE;
+ }
+ else { /* one is LT and the other GT */
+ *e0 = (*e0)->next;
+ *e1 = (*e1)->next;
+ }
+ }
+ }
+ }
+ else {
+ if ((*e0)->sign < (*e1)->sign)
+ *e0 = (*e0)->next;
+ else
+ *e1 = (*e1)->next;
+ }
+ }
+ return FALSE;
+}
+
+/*
+ * Given a vector, this function returns a pointer to the table entry (a list)
+ * which classifies the vector in terms of what signs can be emitted for
+ * that vector.
+ */
+static TblEntry *
+setList(Vector v)
+{
+ int sign_a, sign_b;
+ int tmp;
+
+ sign_a = mpz_sgn(v[0]);
+ sign_b = mpz_sgn(v[1]);
+ switch (sign_a) {
+ case 0 :
+ switch (sign_b) {
+ case 0 : /* a == 0, b == 0 */
+ return &tee1;
+ case 1 : /* a == 0, b > 0 */
+ return &teg1;
+ case -1 : /* a == 0, b < 0 */
+ return &tel1;
+ default :
+ Error(FATAL,E_INT,"signFromVector",
+ "bad value returned by mpz_sgn");
+ return NULL;
+ }
+ break;
+ case 1 :
+ switch (sign_b) {
+ case 0 : /* a > 0, b == 0 */
+ return &tge1;
+
+ case 1 : /* a > 0, b > 0 */
+ tmp = mpz_cmp(v[0], v[1]);
+ switch (MPZ_SIGN(tmp)) {
+ case 0 : /* a > 0, b > 0, a == b */
+ return &tgge1;
+ case 1 : /* a > b, b > 0, a > b */
+ return &tggg1;
+ case -1 : /* a > b, b > 0, a < b */
+ return &tggl1;
+ }
+ break;
+
+ case -1 : /* a > 0, b < 0 */
+ /* it is faster to negate a twice */
+ mpz_neg(v[0], v[0]);
+ tmp = mpz_cmp(v[0], v[1]);
+ mpz_neg(v[0], v[0]);
+ switch (MPZ_SIGN(tmp)) {
+ case 0 : /* a > 0, b < 0, -a == b */
+ return &tgle1;
+ case 1 : /* a > b, b < 0, -a > b */
+ return &tglg1;
+ case -1 : /* a > b, b < 0, -a < b */
+ return &tgll1;
+ }
+ break;
+
+ default :
+ Error(FATAL,E_INT,"signFromVector",
+ "bad value returned by mpz_sgn");
+ return NULL;
+ }
+ break;
+
+ case -1 :
+ switch (sign_b) {
+ case 0 : /* a < 0, b == 0 */
+ return &tle1;
+ case 1 : /* a < 0, b > 0 */
+ /* it is faster to negate a twice */
+ mpz_neg(v[0], v[0]);
+ tmp = mpz_cmp(v[0], v[1]);
+ mpz_neg(v[0], v[0]);
+ switch (MPZ_SIGN(tmp)) {
+ case 0 : /* a < 0, b > 0, -a == b */
+ return &tlge1;
+ case 1 : /* a < b, b > 0, -a > b */
+ return &tlgg1;
+ case -1 : /* a < b, b > 0, -a < b */
+ return &tlgl1;
+ }
+ break;
+
+ case -1 : /* a < 0, b < 0 */
+ tmp = mpz_cmp(v[0], v[1]);
+ switch (MPZ_SIGN(tmp)) {
+ case 0 : /* a < 0, b < 0, a == b */
+ return &tlle1;
+ case 1 : /* a < b, b < 0, a > b */
+ return &tllg1;
+ case -1 : /* a < b, b < 0, a < b */
+ return &tlll1;
+ }
+ break;
+
+ default :
+ Error(FATAL,E_INT,"signFromVector",
+ "bad value returned by mpz_sgn");
+ return NULL;
+ }
+ break;
+
+ default :
+ Error(FATAL,E_INT,"signFromVector","bad value returned by mpz_sgn");
+ return NULL;
+ }
+ return NULL;
+}
+
+static void
+rNoOp(Vector v)
+{
+}
+
+static void
+regSNEG(Vector v)
+{
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+regSZERO(Vector v)
+{
+ mpz_set(v[0], v[1]);
+}
+
+static void
+relSPOS(Vector v)
+{
+ mpz_neg(v[1], v[1]);
+}
+
+static void
+relSNEG(Vector v)
+{
+ mpz_neg(v[1], v[1]);
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+relSZERO(Vector v)
+{
+ mpz_neg(v[1], v[1]);
+ mpz_set(v[0], v[1]);
+}
+
+static void
+rgeSNEG(Vector v)
+{
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+rgeSINF(Vector v)
+{
+ mpz_set(v[1], v[0]);
+}
+
+static void
+rleSPOS(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+}
+
+static void
+rleSNEG(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+rleSINF(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ mpz_set(v[1], v[0]);
+}
+
+static void
+rggeSZERO(Vector v)
+{
+ mpz_mul_2exp(v[0], v[0], (unsigned long) 1);
+ mpz_set_ui(v[1], (unsigned long) 0);
+}
+
+static void
+rggeSINF(Vector v)
+{
+ mpz_mul_2exp(v[1], v[0], (unsigned long) 1);
+ mpz_set_ui(v[0], (unsigned long) 0);
+}
+
+static void
+rgggSINF(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_sub(v[0], v[0], v[1]);
+ mpz_add(v[1], tmp_z, v[1]);
+}
+
+static void
+rgglSZERO(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_add(v[0], v[0], v[1]);
+ mpz_sub(v[1], v[1], tmp_z);
+}
+
+static void
+rll(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ mpz_neg(v[1], v[1]);
+}
+
+static void
+rlleSZERO(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ mpz_mul_2exp(v[0], v[0], (unsigned long) 1);
+ mpz_set_ui(v[1], (unsigned long) 0);
+}
+
+static void
+rlleSINF(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ mpz_mul_2exp(v[1], v[0], (unsigned long) 1);
+ mpz_set_ui(v[0], (unsigned long) 0);
+}
+
+static void
+rllgSZERO(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_add(v[0], v[0], v[1]);
+ mpz_neg(v[0], v[0]);
+ mpz_sub(v[1], tmp_z, v[1]);
+}
+
+static void
+rlllSINF(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_sub(v[0], v[1], v[0]);
+ mpz_add(v[1], tmp_z, v[1]);
+ mpz_neg(v[1], v[1]);
+}
+
+static void
+rgl(Vector v)
+{
+ mpz_neg(v[1], v[1]);
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+rgleSZERO(Vector v)
+{
+ mpz_mul_2exp(v[1], v[0], (unsigned long) 1);
+ mpz_set_ui(v[0], (unsigned long) 0);
+}
+
+static void
+rgleSINF(Vector v)
+{
+ mpz_mul_2exp(v[0], v[0], (unsigned long) 1);
+ mpz_set_ui(v[1], (unsigned long) 0);
+}
+
+static void
+rglgSZERO(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_add(v[0], v[0], v[1]);
+ mpz_neg(v[0], v[0]);
+ mpz_sub(v[1], tmp_z, v[1]); /* - (b - a) */
+}
+
+static void
+rgllSINF(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_sub(v[0], v[0], v[1]);
+ mpz_add(v[1], tmp_z, v[1]);
+}
+
+static void
+rlg(Vector v)
+{
+ mpz_neg(v[0], v[0]);
+ MPZ_SWAP(v[0], v[1]);
+}
+
+static void
+rlgeSZERO(Vector v)
+{
+ mpz_mul_2exp(v[1], v[1], (unsigned long) 1);
+ mpz_set_ui(v[0], (unsigned long) 0);
+}
+
+static void
+rlgeSINF(Vector v)
+{
+ mpz_mul_2exp(v[0], v[1], (unsigned long) 1);
+ mpz_set_ui(v[1], (unsigned long) 0);
+}
+
+static void
+rlggSINF(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_sub(v[0], v[1], v[0]);
+ mpz_add(v[1], tmp_z, v[1]);
+ mpz_neg(v[1], v[1]);
+}
+
+static void
+rlglSZERO(Vector v)
+{
+ mpz_set(tmp_z, v[0]);
+ mpz_add(v[0], v[0], v[1]);
+ mpz_sub(v[1], v[1], tmp_z);
+}
generated by cgit v1.2.3 (git 2.25.1) at 2024-09-19 16:21:56 +0000