implementing some new functionality to bitvecotr

bitvector.cc

@@ -101,9 +101,7 @@ return r;
 
 bitvector& bitvector::operator&=(const bitvector &rhs)
 {
-#ifdef DEBUG
+if(size()<rhs.size()) resize(rhs.size(),true);
-if(nn!=rhs.nn || modulo!=rhs.modulo) laerror("operation on incompatible bitvectors");
-#endif
 copyonwrite();
 for(int i=0; i<nn; ++i) v[i] &= rhs.v[i];
 return *this;
@@ -111,9 +109,7 @@ return *this;
 
 bitvector& bitvector::operator|=(const bitvector &rhs)
 {
-#ifdef DEBUG
+if(size()<rhs.size()) resize(rhs.size(),true);
-if(nn!=rhs.nn || modulo!=rhs.modulo) laerror("operation on incompatible bitvectors");
-#endif
 copyonwrite();
 for(int i=0; i<nn; ++i) v[i] |= rhs.v[i];
 return *this;
@@ -121,9 +117,7 @@ return *this;
 
 bitvector& bitvector::operator^=(const bitvector &rhs)
 {
-#ifdef DEBUG
+if(size()<rhs.size()) resize(rhs.size(),true);
-if(nn!=rhs.nn || modulo!=rhs.modulo) laerror("operation on incompatible bitvectors");
-#endif
 copyonwrite();
 for(int i=0; i<nn; ++i) v[i] ^= rhs.v[i];
 return *this;
@@ -157,6 +151,61 @@ return s;
 #endif
 
 
+bitvector& bitvector::operator>>=(unsigned int i)
+{
+if(i==0) return *this;
+copyonwrite();
+unsigned int imod = i%blockbits;
+unsigned int ishift = i/blockbits;
+for(int dest=0; dest<nn; ++dest)
+	{
+	int src=dest+ishift;
+	if(src>=nn) v[dest]=0;
+	else
+		{
+		v[dest] = v[src]>>imod;
+		if(imod && (src+1<nn)) v[dest] |= (v[src+1]&((1ULL<<imod)-1)) <<(blockbits-imod);
+		}
+	}
+return *this;
+}
+
+bitvector& bitvector::leftshift(unsigned int i, bool autoresize)
+{
+if(i==0) return *this;
+copyonwrite();
+unsigned int imod = i%blockbits;
+unsigned int ishift = i/blockbits;
+if(autoresize) resize(size()+i,true);
+for(int dest=nn-1; dest>=0; --dest)
+        {
+        int src=dest-ishift;
+        if(src<0) v[dest]=0;
+        else
+                {
+                v[dest] = v[src]<<imod;
+                if(imod && (src-1>=0)) v[dest] |= (v[src-1]& (((1ULL<<imod)-1) <<(blockbits-imod)))>>(blockbits-imod);
+                }
+        }
+
+return *this;
+}
+
+
+void  bitvector::randomize()
+{
+copyonwrite();
+for(int i=0; i<nn; ++i) v[i]=RANDINT64();
+//zero the excess bits
+if(modulo) 
+	{
+	bitvector_block mask = (1ULL<<modulo)-1;
+	v[nn-1] &= mask;
+	}
+}
+
+
+
 unsigned int bitvector::population(const unsigned int before) const
 {
 if(before) laerror("before parameter in population() not implemented yet");

bitvector.h

@@ -20,13 +20,14 @@
 #define _BITVECTOR_H_
 
 #include "vec.h"
+#include <stdint.h>
+
 
 namespace LA {
-
 //compressed storage of large bit vectors
-//any reasonable compiler changes the dividions and modulos to shift/and instructions
+//let's now use 64-bit blocks exclusively 
 
-typedef unsigned long bitvector_block; //should be automatically portable and efficiently use wordlength of each machine (32 vs 64)
+typedef uint64_t bitvector_block; 
 
 #define blockbits (8*sizeof(bitvector_block))
 
@@ -47,19 +48,21 @@ public:
 	//operator= seems to be correctly synthetized by the compiler
 	//override dereferencing to address single bits, is however possible
 	//only in the const context (otherwise we would have to define a type which, when assigned to, changes a single bit - possible but probably inefficient)
-	void resize(const unsigned int n) {NRVec<bitvector_block>::resize((n+blockbits-1)/blockbits); modulo=n%blockbits;};
+	void resize(const unsigned int n, bool preserve=false) {NRVec<bitvector_block>::resize((n+blockbits-1)/blockbits,preserve); modulo=n%blockbits;};
 	unsigned int size() const {return (nn*blockbits)-blockbits+(modulo?modulo:blockbits);};
 	//arguments must be unsigned to keep the resulting assembly code simple and efficient
         const bool operator[](const unsigned int i) const {return (v[i/blockbits] >>(i%blockbits))&1UL;};
+	const bool get(const unsigned int i) const {return (*this)[i];};
 	bitvector_block getblock(const unsigned int i) const {return v[i];}; //integer interpretation
 	void setblock(const unsigned int i, const bitvector_block b) {v[i]=b;};
 	int getblocksize() const {return 8*sizeof(bitvector_block);};
 	void set(const unsigned int i) {v[i/blockbits] |= (1UL<<(i%blockbits));};
 	void reset(const unsigned int i) {v[i/blockbits] &= ~(1UL<<(i%blockbits));};
-	const bool get(const unsigned int i) {return (v[i/blockbits] >>(i%blockbits))&1UL;};
 	const bool assign(const unsigned int i, const bool r) {if(r) set(i); else reset(i); return r;};
 	void clear() {copyonwrite(true); memset(v,0,nn*sizeof(bitvector_block));};
 	void fill() {memset(v,0xff,nn*sizeof(bitvector_block));};
+	bool iszero() const {for(int i=0; i<nn; ++i) if(v[i]) return false; return true;};
+	void randomize();
 	bool operator!=(const bitvector &rhs) const;
 	bool operator==(const bitvector &rhs) const {return !(*this != rhs);};
 	bool operator>(const bitvector &rhs) const;
@@ -70,13 +73,22 @@ public:
 	bitvector& operator&=(const bitvector &rhs);
 	bitvector& operator|=(const bitvector &rhs);
 	bitvector& operator^=(const bitvector &rhs);
+	bitvector& operator+=(const bitvector &rhs) {return (*this)^=rhs;}; //addition modulo 2
+	bitvector& operator-=(const bitvector &rhs) {return (*this)^=rhs;}; //subtraction modulo 2
 	bitvector operator&(const bitvector &rhs) const {return bitvector(*this) &= rhs;};
 	bitvector operator|(const bitvector &rhs) const {return bitvector(*this) |= rhs;};
 	bitvector operator^(const bitvector &rhs) const {return bitvector(*this) ^= rhs;};
+	bitvector operator+(const bitvector &rhs) const {return *this ^ rhs;}; //addition modulo 2
+	bitvector operator-(const bitvector &rhs) const {return *this ^ rhs;}; //subtraction modulo 2
         unsigned int operator%(const bitvector &y) const; //number of differing bits
 	unsigned int population(const unsigned int before=0) const; //number of 1's
 	//extended, truncated const i.e. not on *this but return new entity, take care of modulo's bits
-	//logical shifts  <<= >>= << >> not implemented yet
+	//logical shifts
+	bitvector& operator>>=(unsigned int i);
+	bitvector& leftshift(unsigned int i, bool autoresize=false);
+	bitvector& operator<<=(unsigned int i) {return leftshift(i,true);};
+	bitvector operator>>(unsigned int i) {bitvector r(*this); return r>>=i;};
+	bitvector operator<<(unsigned int i) {bitvector r(*this); return r<<=i;};
 	//logical rotations not implemented yet
 	//unformatted file IO
 	void read(int fd, bool dimensions=1, bool transp=0);

la_random.cc

@@ -6,5 +6,6 @@ namespace LA {
 WEAK_SYMBOL  double randdouble() {return random()/(1.+RAND_MAX);}
 WEAK_SYMBOL  double randdoublesigned() {return 2.*random()/(1.+RAND_MAX)-1.;}
 WEAK_SYMBOL  int randint32() {return random();}
+WEAK_SYMBOL  uint64_t randint64() {uint64_t r = random(); r<<=32; r|= random(); return r; }
 }//namespace
 

la_random.h

@@ -1,11 +1,14 @@
 #ifndef _LA_RANDOM_H
 #define _LA_RANDOM_H
 
+#include <stdint.h>
+
 namespace LA {
 
 extern double randdouble();
 extern double randdoublesigned();
 extern int randint32();
+extern uint64_t randint64();
 
 //RANDOM numbers defaulting to standard library but switchable to user's functions
 #ifndef RANDDOUBLE
@@ -20,6 +23,11 @@ extern int randint32();
 #define RANDINT32 LA::randint32
 #endif
 
+#ifndef RANDINT64
+#define RANDINT64 LA::randint64
+#endif
+
+
 #ifdef __GNUC__
 #define WEAK_SYMBOL __attribute__((weak)) 
 #else

t.cc

@@ -2784,7 +2784,7 @@ NRSMat<char> adjperm = adj.permuted(p);
 cout <<"resorted graph = "<<adjperm<<endl;
 }
 
-if(1)
+if(0)
 {
 int seed;
 int f=open("/dev/random",O_RDONLY);
@@ -2852,4 +2852,29 @@ cout <<endl<<"Inverse via svd\n"<<ainv2<<endl;
 cout <<"Difference of inverses = "<<(ainv-ainv2).norm()<<endl;
 }
 
+if(1)
+{
+int seed;
+int f=open("/dev/random",O_RDONLY);
+if(sizeof(int)!=read(f,&seed,sizeof(int))) laerror("cannot read /dev/random");
+close(f);
+srand(seed);
+
+int n;
+cin >>n;
+bitvector v(n);
+v.randomize();
+//do{
+//	cout <<v <<endl;
+//	v>>=1;
+//}while(!v.iszero());
+
+
+for(int i=0; i<n; ++i) 
+	{
+	cout <<v <<endl;
+        v<<=1;
+	}
+}
+
 }