First cut at starUpperBound algorithm.

Author: highperformancecoder

classdesc

@@ -127,8 +127,9 @@ ifdef AEGIS
   # for now, we maintain the code for C++98 standard
   # noexcept-type is a warning about ABI changes coming in C++17
   # Travis complains about unused-local-typedefs, so shut it up.
-  CXXFLAGS+=-Wno-noexcept-type -Wno-unused-local-typedefs -std=c++20
+  CXXFLAGS+=-Wno-noexcept-type -Wno-unused-local-typedefs
 endif
+CXXFLAGS+=-std=c++20
 
 ifndef MXE
   # needed for generating .so files

include/Makefile

@@ -1,5 +1,7 @@
 #include "starComplexity.h"
+#ifdef SYCL_LANGUAGE_VERSION
 #include "vecBitSet.cd"
+#endif
 #include "starComplexity.cd"
 #include "netcomplexity.h"
 #include "pythonBuffer.h"
@@ -33,6 +35,8 @@ linkRep edge(unsigned i, unsigned j)
   return linkRep(1)<<(i*(i-1)/2+j);
 }
 
+
+
 void StarComplexityGen::generateElementaryStars(unsigned nodes)
 {
   for (unsigned i=0; i<nodes; ++i)
@@ -431,7 +435,7 @@ NautyRep toNautyRep(linkRep g, unsigned nodes)
   NautyRep n(nodes);
   for (unsigned i=0; i<nodes; ++i)
     for (unsigned j=0; j<i; ++j)
-      if (g&(linkRep(1)<<(i*(i-1)/2+j)))
+      if (!(g&(linkRep(1)<<(i*(i-1)/2+j))).empty())
         n(i,j)=n(j,i)=true;
   return n;
 }
@@ -545,3 +549,59 @@ GraphComplexity StarComplexityGen::complexity(linkRep g) const
     }
   return r;
 }
+
+
+unsigned starUpperBound(linkRep x, unsigned nodes) 
+{
+  // Firstly remove those nodes that are full stars
+  vector<unsigned> fullStars;
+  for (unsigned i=0; i<nodes; ++i)
+    {
+      bool fullStar=true;
+      for (unsigned j=0; fullStar && j<nodes; ++j)
+        if (i!=j && !x(i,j))
+          fullStar=false;
+      if (fullStar) fullStars.push_back(i);
+    }
+
+  unsigned stars=fullStars.size();
+  for (auto i: fullStars)
+    for (unsigned j=0; j<nodes; ++j)
+      x&=~edge(i,j); // remove all edges to this node
+
+
+  // calculate node degree
+  map<unsigned,unsigned> nodeDegree;
+  for (unsigned i=0; i<maxNodes; ++i)
+    for (unsigned j=0; j<i; ++j)
+      if (x(i,j))
+        {
+          ++nodeDegree[i];
+          ++nodeDegree[j];
+        }
+
+
+  struct SecondLess
+  {
+    bool operator()(const pair<unsigned,unsigned>& x, const pair<unsigned,unsigned>& y)
+      const {return x.second < y.second;}
+  };
+  
+  // compute the number of operations xᵢ∪(xₖ∩xₗ…)
+  auto maxDegree=max_element(nodeDegree.begin(), nodeDegree.end(), SecondLess());
+
+  while (maxDegree->second>0)
+    {
+      stars+=1+maxDegree->second;
+      maxDegree->second=0; // accounted for all links to this node
+      for (unsigned j=0; j<maxNodes; ++j)
+        if (x(maxDegree->first, j) && nodeDegree[j])
+          --nodeDegree[j];
+      maxDegree=max_element(nodeDegree.begin(), nodeDegree.end(), SecondLess());
+    }
+
+  return stars;
+}
+
+unsigned StarComplexityGen::starUpperBound(const linkRep& x) const
+{return min(::starUpperBound(x, elemStars.size()), ::starUpperBound(~x,elemStars.size()));}

models/starComplexity.cc

@@ -1,5 +1,5 @@
 // hard code maximum number of nodes
-constexpr unsigned maxNodes=22, maxStars=2*maxNodes-1;
+constexpr unsigned maxNodes=8, maxStars=2*maxNodes-1;
 
 //using linkRep=unsigned long long;
 
@@ -19,7 +19,7 @@ class linkRepImpl
   CLASSDESC_ACCESS(linkRepImpl);
 public:
   linkRepImpl()=default;
-  linkRepImpl(unsigned long long x) {
+  linkRepImpl(unsigned x) {
     static_assert(sizeof(data)>=sizeof(x));
     memcpy(data,&x,sizeof(x));
     memset(((char*)data)+sizeof(x),0,sizeof(data)-sizeof(x));
@@ -47,12 +47,22 @@ class linkRepImpl
       r.data[i+d.quot]=data[i]<<d.rem;
     return r;
   }
-  operator const void*() const {
+  // return true if there is a link between nodes i and j
+  bool operator()(unsigned i,unsigned j) const {
+        if (i<j) std::swap(i,j);
+        auto d=div(i*(i-1)/2+j, int(8*sizeof(Impl)));
+        return (data[d.quot] & (Impl(1)<<d.rem))!=0;
+  }
+  bool empty() const {
     for (unsigned i=0; i<size; ++i)
-      if (data[i]) return data;
-    return nullptr;
+      if (data[i]) return false;
+    return true;
+  }
+  bool operator==(const linkRepImpl& x) const {
+    for (unsigned i=0; i<size; ++i)
+      if (data[i]!=x.data[i]) return false;
+    return true;
   }
-
   auto operator<=>(const linkRepImpl& x) const {
     for (unsigned i=0; i<size; ++i)
       if (auto r=data[i]<=>x.data[i]; r!=0)
@@ -134,5 +144,8 @@ struct StarComplexityGen
   linkRep complement(linkRep) const;
   unsigned symmStar(linkRep) const;
   GraphComplexity complexity(linkRep) const;
+
+  /// return an upper bound on the number of stars in the link representation
+  unsigned starUpperBound(const linkRep&) const;
 };
 

models/starComplexity.h

@@ -3,7 +3,7 @@
 
 from starComplexity import starC
 from datetime import datetime
-nodes=22
+nodes=7
 # computed from max_{l\in[0,L]}min(n+L-l,2l) where L=n(n-1)/2
 maxStars=0
 L=int(0.5*nodes*(nodes-1))
@@ -12,7 +12,7 @@
     maxStars=max(maxStars,v)
 
 print('maxStars=',maxStars)
-maxStars=12
+maxStars=6
 
 #starC.blockSize(256)
 starC.blockSize(4096)
@@ -21,7 +21,7 @@
 #starC.blockSize(40320)
 
 starC.generateElementaryStars(nodes)
-for numStars in range(2,maxStars+1):
+for numStars in range(1,maxStars+1):
     starC.fillStarMap(numStars)
     print('completed',numStars,datetime.now())
     starC.canonicaliseStarMap()
@@ -33,5 +33,7 @@
             links=0
             for j in i:
                 links+=bin(j).count('1')
-            print(id,links,starC.symmStar(i)-1,c.complexity(),c.starComplexity(),starC.counts[i],starC.counts[starC.complement(i)()],sep=',',file=out)
+
+            #print(id,links,starC.symmStar(i)-1,c.complexity(),c.starComplexity(),starC.counts[i],sep=',',file=out)
+            print(id,bin(i[0]),starC.starMap[i],starC.starUpperBound(i));
             id+=1