[mathicgb] 214/393: Added functions for testing a*b==c to MonoMonoid and forwarding from PolyRing.

Doug Torrance dtorrance-guest at moszumanska.debian.org
Fri Apr 3 15:59:03 UTC 2015


This is an automated email from the git hooks/post-receive script.

dtorrance-guest pushed a commit to branch upstream
in repository mathicgb.

commit 15092b8ccd356a120efbd059a31951feb0060715
Author: Bjarke Hammersholt Roune <bjarkehr.code at gmail.com>
Date:   Sun Mar 31 00:11:15 2013 +0100

    Added functions for testing a*b==c to MonoMonoid and forwarding from PolyRing.
---
 src/mathicgb/MonoMonoid.hpp | 100 +++++++++++++++++++++++++++++++++++++++++++-
 src/mathicgb/PolyRing.hpp   |  20 +++++++++
 src/test/MonoMonoid.cpp     |  24 +++++++++--
 3 files changed, 140 insertions(+), 4 deletions(-)

diff --git a/src/mathicgb/MonoMonoid.hpp b/src/mathicgb/MonoMonoid.hpp
index 9e468ec..9e7a070 100755
--- a/src/mathicgb/MonoMonoid.hpp
+++ b/src/mathicgb/MonoMonoid.hpp
@@ -8,6 +8,7 @@
 #include <istream>
 #include <ostream>
 #include <cstdlib>
+#include <cstring>
 #include <mathic.h>
 
 /// Implements the monoid of (monic) monomials with integer
@@ -187,8 +188,103 @@ public:
     return static_cast<HashValue>(access(mono, hashIndex()));
   }
 
+  /// Returns true if a and b are equal. Includes check for component.
   bool equal(ConstMonoRef a, ConstMonoRef b) const {
-    return std::equal(begin(a), end(a), begin(b));
+    for (auto i = componentIndex(); i != exponentsIndexEnd(); ++i)
+      if (access(a, i) != access(b, i))
+        return false;
+    return true;
+  }
+
+  /// As equal(), but optimized for the case where true is returned.
+  bool equalHintTrue(ConstMonoRef a, ConstMonoRef b) const {
+    // if a[i] != b[i] then a[i] ^ b[i] != 0, so the or of all xors is zero
+    // if and only if a equals b. This way we avoid having a branch to check
+    // equality for every iteration of the loop, which is a win in the case
+    // that none of the early-exit branches are taken - that is, when a equals
+    // b.
+    Exponent orOfXor = 0;
+    for (VarIndex i = lastExponentIndex(); i != beforeComponentIndex(); --i)
+      orOfXor |= access(a, i) ^ access(b, i);
+    MATHICGB_ASSERT((orOfXor == 0) == equal(a, b));
+    return orOfXor == 0;
+  }
+
+  bool isProductOf(
+    ConstMonoRef a,
+    ConstMonoRef b,
+    ConstMonoRef ab
+  ) const {
+    for (VarIndex i = componentIndex(); i != exponentsIndexEnd(); ++i)
+      if (access(ab, i) != access(a, i) + access(b, i))
+        return false;
+    return true;
+  }
+
+  bool isProductOfHintTrue(
+    ConstMonoRef a, 
+    ConstMonoRef b, 
+    ConstMonoRef ab
+  ) const {
+    // We compare more than one exponent at a time using 64 bit integers. This 
+    // might go one 32 bit value at the end too far, but since that space is
+    // either a degree or a hash value that is fine --- those values will also
+    // match if the monomials are equal. This does not work for negative
+    // exponents since the overflowing bit will go into the next word.
+    // It is OK that the degree field can be negative (a field we might go
+    // into without caring about it because it shares a 64 bit field with
+    // the last exponent), because it is at the end so the overflowing
+    // bit will not interfere.
+    
+    // todo: ensure 8 byte alignment. Though there seem to be no ill effects
+    // for unaligned access. Performance seems to be no worse than for using
+    // 32 bit integers directly.
+
+    if (sizeof(Exponent) != 4)
+      return isProductOf(a, b, ab);
+
+    uint64 orOfXor = 0;
+    for (VarIndex i = varCount() / 2; i != beforeComponentIndex(); --i) {
+      MATHICGB_ASSERT(access(a, i*2) >= 0);
+      MATHICGB_ASSERT(i == varCount() / 2 || access(a, i*2+1) >= 0);
+      
+      uint64 A, B, AB;
+      // We have to use std::memcpy here because just casting to a int64 breaks
+      // the strict aliasing rule which implies undefined behavior. Both MSVC and
+      // gcc don't actually call memcpy here. MSVC is a tiny bit slower for this
+      // code than for casting while GCC seems to be exactly the same speed.
+      std::memcpy(&A, ptr(a, i*2), 8);
+      std::memcpy(&B, ptr(b, i*2), 8);
+      std::memcpy(&AB, ptr(ab, i*2), 8);
+      orOfXor |= AB ^ (A + B);
+    }
+    MATHICGB_ASSERT((orOfXor == 0) == isProductOf(a, b, ab));
+    return orOfXor == 0; 
+  }
+
+  MATHICGB_INLINE bool isTwoProductsOfHintTrue(
+    ConstMonoRef a1,
+    ConstMonoRef a2,
+    ConstMonoRef b,
+    ConstMonoRef a1b,
+    ConstMonoRef a2b
+  ) const {
+    if (sizeof(Exponent) != 4)
+      return (isProductOf(a1, b, a1b) && isProductOf(a2, b, a2b));
+
+    uint64 orOfXor = 0;
+    for (VarIndex i = varCount() / 2; i != beforeComponentIndex(); --i) {
+      uint64 A1, A2, B, A1B, A2B;
+      std::memcpy(&A1, ptr(a1, i*2), 8);
+      std::memcpy(&A2, ptr(a2, i*2), 8);
+      std::memcpy(&B, ptr(b, i*2), 8);
+      std::memcpy(&A1B, ptr(a1b, i*2), 8);
+      std::memcpy(&A2B, ptr(a2b, i*2), 8);
+      orOfXor |= (A1B ^ (A1 + B)) | (A2B ^ (A2 + B));
+    }
+    MATHICGB_ASSERT
+      ((orOfXor == 0) == (isProductOf(a1, b, a1b) && isProductOf(a2, b, a2b)));
+    return orOfXor == 0;
   }
 
   /// Returns the hash of the product of a and b.
@@ -931,6 +1027,8 @@ private:
 
   VarIndex lastEntryIndex() const {return hashIndex();}
   VarIndex beforeFirstIndex() const {return static_cast<VarIndex>(-1);}
+  VarIndex lastExponentIndex() const {return orderIndexEnd() - 1;}
+  VarIndex beforeComponentIndex() const {return componentIndex() - 1;}
 
   const VarIndex mVarCount;
   const VarIndex mOrderEntryCount;
diff --git a/src/mathicgb/PolyRing.hpp b/src/mathicgb/PolyRing.hpp
index bc176f7..6ff481c 100755
--- a/src/mathicgb/PolyRing.hpp
+++ b/src/mathicgb/PolyRing.hpp
@@ -583,16 +583,23 @@ inline exponent PolyRing::weight(ConstMonomial a) const {
 
 inline bool PolyRing::monomialEQ(ConstMonomial a, ConstMonomial b) const
 {
+#ifdef MATHICGB_USE_MONOID
+  return monoid().equal(a, b);
+#else
   for (size_t i = 0; i <= mNumVars; ++i)
     if (a[i] != b[i])
       return false;
   return true;
+#endif
 }
 
 inline bool PolyRing::monomialEqualHintTrue(
   const ConstMonomial a,
   const ConstMonomial b
 ) const {
+#ifdef MATHICGB_USE_MONOID
+  return monoid().equalHintTrue(a, b);
+#else
   // if a[i] != b[i] then a[i] ^ b[i] != 0, so the or of all xors is zero
   // if and only if a equals b. This way we avoid having a branch to check
   // equality for every iteration of the loop, which is a win in the case
@@ -603,6 +610,7 @@ inline bool PolyRing::monomialEqualHintTrue(
   const bool areEqual = (orOfXor == 0);
   MATHICGB_ASSERT(areEqual == monomialEQ(a, b));
   return areEqual;
+#endif
 }
 
 inline bool PolyRing::monomialIsProductOfHintTrue(
@@ -610,6 +618,9 @@ inline bool PolyRing::monomialIsProductOfHintTrue(
   const ConstMonomial b, 
   const ConstMonomial ab
 ) const {
+#ifdef MATHICGB_USE_MONOID
+  return monoid().isProductOfHintTrue(a, b, ab);
+#else
   // We compare more than one exponent at a time using 64 bit integers. This 
   // might go one 32 bit value at the end too far, but since that space is
   // either a degree or a hash value that is fine --- those values will also
@@ -645,6 +656,7 @@ inline bool PolyRing::monomialIsProductOfHintTrue(
   MATHICGB_ASSERT((orOfXor == 0) == monomialIsProductOf(a, b, ab));
 
   return orOfXor == 0; 
+#endif
 }
 
 MATHICGB_INLINE bool PolyRing::monomialIsTwoProductsOfHintTrue(
@@ -654,6 +666,9 @@ MATHICGB_INLINE bool PolyRing::monomialIsTwoProductsOfHintTrue(
   const ConstMonomial a1b,
   const ConstMonomial a2b
 ) const {
+#ifdef MATHICGB_USE_MONOID
+  return monoid().isTwoProductsOfHintTrue(a1, a2, b, a1b, a2b);
+#else
   if (sizeof(exponent) < 4)
     return (monomialIsProductOf(a1, b, a1b) &&
       monomialIsProductOf(a2, b, a2b));
@@ -672,6 +687,7 @@ MATHICGB_INLINE bool PolyRing::monomialIsTwoProductsOfHintTrue(
     (monomialIsProductOf(a1, b, a1b) && monomialIsProductOf(a2, b, a2b)));
 
   return orOfXor == 0;
+#endif
 }
 
 inline bool PolyRing::monomialIsProductOf(
@@ -679,10 +695,14 @@ inline bool PolyRing::monomialIsProductOf(
   ConstMonomial b, 
   ConstMonomial ab
 ) const {
+#ifdef MATHICGB_USE_MONOID
+  return monoid().isProductOf(a, b, ab);
+#else
   for (size_t i = 0; i <= mNumVars; ++i)
     if (ab[i] != a[i] + b[i])
       return false;
   return true;
+#endif
 }
 
 inline void PolyRing::monomialMult(ConstMonomial a, 
diff --git a/src/test/MonoMonoid.cpp b/src/test/MonoMonoid.cpp
index 241fb5d..3d036fd 100755
--- a/src/test/MonoMonoid.cpp
+++ b/src/test/MonoMonoid.cpp
@@ -322,6 +322,12 @@ TEST(MonoMonoid, MultiplyDivide) {
     ASSERT_EQ(m.hashOfProduct(a, b), m.hash(c));
     ASSERT_EQ(m.hashOfProduct(a, b), m.hashOfProduct(b, a));
 
+    // isProductOf
+    ASSERT_TRUE(m.isProductOf(a, b, c));
+    ASSERT_TRUE(m.isProductOfHintTrue(a, b, c));
+    ASSERT_TRUE(m.isTwoProductsOfHintTrue(a, a, b, c, c));
+    
+
     // a*b == c using multiply
     m.multiply(a, b, mono);
     ASSERT_TRUE(m.equal(c, mono));
@@ -359,6 +365,12 @@ TEST(MonoMonoid, MultiplyDivide) {
       ASSERT_FALSE(m.lessThan(mono, b));
       ASSERT_TRUE(m.compare(mono, b) == Monoid::GreaterThan);
       ASSERT_FALSE(m.divides(mono, b));
+
+      ASSERT_FALSE(m.isProductOf(a, c, b));
+      ASSERT_FALSE(m.isProductOfHintTrue(a, c, b));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(c, c, a, b, b));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(b, c, a, c, b));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(c, b, a, b, c));
     } else {
       ASSERT_TRUE(m.equal(b, mono));
       ASSERT_TRUE(m.compare(b, mono) == Monoid::EqualTo);
@@ -370,6 +382,12 @@ TEST(MonoMonoid, MultiplyDivide) {
       ASSERT_FALSE(m.lessThan(mono, a));
       ASSERT_TRUE(m.compare(mono, a) == Monoid::GreaterThan);
       ASSERT_FALSE(m.divides(mono, a));
+
+      ASSERT_FALSE(m.isProductOf(c, b, a));
+      ASSERT_FALSE(m.isProductOfHintTrue(b, c, a));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(c, c, b, a, a));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(a, c, b, c, a));
+      ASSERT_FALSE(m.isTwoProductsOfHintTrue(c, a, b, a, c));
     } else {
       ASSERT_TRUE(m.equal(a, mono));
       ASSERT_TRUE(m.compare(a, mono) == Monoid::EqualTo);
@@ -391,13 +409,13 @@ TEST(MonoMonoid, MultiplyDivide) {
     ASSERT_TRUE(m.equal(mono, a));
   };
   check("1 1 1");
-  check("a 1 a");
+  check("a<5> 1 a<5>");
   check("1 Vx Vx");
   check("aV bx abxV");
   check("a a2 a3");
-  check("V V2 V3");
+  check("V<2> V2 V3<2>");
   check("arlgh svug arlg2hsvu");
-  check("abcdefghiV ab2c3d4e5f6g7h8i9V11 a2b3c4d5e6f7g8h9i10V12");
+  check("abcdefghiV<7> ab2c3d4e5f6g7h8i9V11 a2b3c4d5e6f7g8h9i10V12<7>");
 }
 
 TEST(MonoMonoid, Order) {

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/mathicgb.git



More information about the debian-science-commits mailing list