From d4c7c27f5dea5fca5fbb660a7931773c5bbf2cfe Mon Sep 17 00:00:00 2001 From: Meghana Sistla Date: Sat, 18 Oct 2025 17:29:08 +0000 Subject: [PATCH 1/3] Made changes for equivalence checking problem --- CFLOBDD/Makefile | 68 ++ CFLOBDD/main.cpp | 3 +- CFLOBDD/matmult_map.cpp | 2 +- CFLOBDD/matrix1234_complex_float_boost.cpp | 28 + CFLOBDD/matrix1234_complex_float_boost.h | 12 + ...atrix1234_complex_float_boost_top_node.cpp | 141 +++- .../matrix1234_complex_float_boost_top_node.h | 101 +++ CFLOBDD/matrix1234_node.cpp | 33 +- CFLOBDD/tests_cfl.cpp | 177 +++++ CFLOBDD/tests_cfl.h | 2 + CFLOBDD/vector_complex_float_boost.cpp | 7 + CFLOBDD/vector_complex_float_boost.h | 1 + .../vector_complex_float_boost_top_node.cpp | 14 + CFLOBDD/vector_complex_float_boost_top_node.h | 1 + CFLOBDD/weighted_cflobdd_node_t.cpp | 6 +- CFLOBDD/weighted_matmult_map.cpp | 5 +- CFLOBDD/wmatrix1234_complex_fb_mul.cpp | 20 + CFLOBDD/wmatrix1234_complex_fb_mul.h | 7 + CFLOBDD/wmatrix1234_complex_fb_mul_node.cpp | 721 ++++++++++++------ CFLOBDD/wmatrix1234_complex_fb_mul_node.h | 86 +++ .../wmatrix1234_top_node_complex_fb_mul.cpp | 73 +- CFLOBDD/wmatrix1234_top_node_complex_fb_mul.h | 4 + CFLOBDD/wvector_complex_fb_mul.cpp | 8 + CFLOBDD/wvector_complex_fb_mul.h | 2 + CFLOBDD/wvector_complex_fb_mul_node.cpp | 1 + CFLOBDD/wvector_top_node_complex_fb_mul.cpp | 12 + CFLOBDD/wvector_top_node_complex_fb_mul.h | 1 + 27 files changed, 1285 insertions(+), 251 deletions(-) create mode 100644 CFLOBDD/Makefile diff --git a/CFLOBDD/Makefile b/CFLOBDD/Makefile new file mode 100644 index 0000000..56e3c0a --- /dev/null +++ b/CFLOBDD/Makefile @@ -0,0 +1,68 @@ +# Project Name (executable) +PROJECT = cflobdd +# Compiler +CC = g++-11 + +# Run Options +COMMANDLINE_OPTIONS = #/dev/ttyS0 + +# Compiler options during compilation +COMPILE_OPTIONS = -g -std=c++2a -w -shared -Wall -Wextra -DHAVE_CONFIG_H -Werror -Wunused-but-set-variable +# -ansi -pedantic -Wall + +#Header include directories +HEADERS = -I. -I../../../../boost_1_81_0/ -I. -I.Solver/uwr/bit_vector/ -I.Solver/uwr/assert/ -I.Solver/uwr/matrix/ -I.Solver/uwr/parsing/ +#Libraries for linking +LIBS = + +# Dependency options +DEPENDENCY_OPTIONS = -MM + + +# Subdirs to search for additional source files +SOURCE_FILES := $(shell ls *.cpp) +SOURCE_FILES += $(shell ls Solver/uwr/bit_vector/*.cpp) +SOURCE_FILES += $(shell ls Solver/uwr/parsing/*.cpp) + +# Create an object file of every cpp file +OBJECTS = $(patsubst %.cpp, %.o, $(SOURCE_FILES)) + +# Dependencies +DEPENDENCIES = $(patsubst %.cpp, %.d, $(SOURCE_FILES)) + +dep: + cd cflobdd/cudd-complex-big && make all && cd ../../ + +# Make $(PROJECT) the default target +all: $(PROJECT) +#$(DEPENDENCIES) -shared +$(PROJECT): $(OBJECTS) + $(CC) -o $(PROJECT) $(OBJECTS) $(LIBS) + +# Include dependencies (if there are any) +# ifneq "$(strip $(DEPENDENCIES))" "" +# include $(DEPENDENCIES) +# endif + +# Compile every cpp file to an object +# %.cpp +%.o: %.cpp + $(CC) -c $(COMPILE_OPTIONS) -o $@ $^ $(HEADERS) + +# Build & Run Project +run: $(PROJECT) + ./$(PROJECT) $(COMMANDLINE_OPTIONS) + +# Clean & Debug +.PHONY: makefile-debug +makefile-debug: + +.PHONY: clean +clean: + rm -f $(PROJECT) $(OBJECTS) + +.PHONY: depclean +depclean: + rm -f $(DEPENDENCIES) + +clean-all: clean depclean diff --git a/CFLOBDD/main.cpp b/CFLOBDD/main.cpp index 477bc5e..2f9aad0 100644 --- a/CFLOBDD/main.cpp +++ b/CFLOBDD/main.cpp @@ -40,7 +40,8 @@ int main(int argc, char * argv[]) // Supply a default argument for when invoking from Windows (e.g., for debugging) if (argc == 1) { - std::string default_string = "And"; + // std::string default_string = "And"; + std::string default_string = "testTmp"; CFL_OBDD::CFLTests::runTests(default_string.c_str()); } else { diff --git a/CFLOBDD/matmult_map.cpp b/CFLOBDD/matmult_map.cpp index 6ae81c5..98e74bd 100644 --- a/CFLOBDD/matmult_map.cpp +++ b/CFLOBDD/matmult_map.cpp @@ -37,7 +37,7 @@ // Constructor MatMultMapBody::MatMultMapBody() - : refCount(0), isCanonical(false) + : refCount(0), isCanonical(false), contains_zero_val(false) { hashCheck = NULL; } diff --git a/CFLOBDD/matrix1234_complex_float_boost.cpp b/CFLOBDD/matrix1234_complex_float_boost.cpp index e6cd3fd..0b52379 100644 --- a/CFLOBDD/matrix1234_complex_float_boost.cpp +++ b/CFLOBDD/matrix1234_complex_float_boost.cpp @@ -91,6 +91,24 @@ namespace CFL_OBDD { return CFLOBDD_COMPLEX_BIG(MkSGateInterleavedTop(i)); } + CFLOBDD_COMPLEX_BIG MkSdgGateInterleaved(unsigned int i) + { + assert(i == 1); + return CFLOBDD_COMPLEX_BIG(MkSdgGateInterleavedTop(i)); + } + + CFLOBDD_COMPLEX_BIG MkTGateInterleaved(unsigned int i) + { + assert(i == 1); + return CFLOBDD_COMPLEX_BIG(MkTGateInterleavedTop(i)); + } + + CFLOBDD_COMPLEX_BIG MkTdgGateInterleaved(unsigned int i) + { + assert(i == 1); + return CFLOBDD_COMPLEX_BIG(MkTdgGateInterleavedTop(i)); + } + CFLOBDD_COMPLEX_BIG MkPhaseShiftGateInterleaved(unsigned int i, double theta) { assert(i == 1); @@ -480,6 +498,16 @@ namespace CFL_OBDD { return CFLOBDD_COMPLEX_BIG(MatrixShiftToBConnectionTop(c.root)); } + CFLOBDD_COMPLEX_BIG ConjugateTranspose(CFLOBDD_COMPLEX_BIG c) + { + return CFLOBDD_COMPLEX_BIG(ConjugateTransposeTop(c.root)); + } + + CFLOBDD_COMPLEX_BIG MatrixTranspose(CFLOBDD_COMPLEX_BIG c) + { + return CFLOBDD_COMPLEX_BIG(MatrixTransposeTop(c.root)); + } + // Return the Kronecker product of two matrices CFLOBDD_COMPLEX_BIG KroneckerProduct2Vocs(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2) { diff --git a/CFLOBDD/matrix1234_complex_float_boost.h b/CFLOBDD/matrix1234_complex_float_boost.h index 90c66c7..e37653d 100644 --- a/CFLOBDD/matrix1234_complex_float_boost.h +++ b/CFLOBDD/matrix1234_complex_float_boost.h @@ -43,7 +43,14 @@ namespace CFL_OBDD { typedef CFLOBDD_T CFLOBDD_COMPLEX_BIG; typedef CFLOBDD_T CFLOBDD_FOURIER; + static constexpr double SQRT2_2 = static_cast(0.707106781186547524400844362104849039284835937688474036588L); + static constexpr double PI = static_cast(3.141592653589793238462643383279502884197169399375105820974L); + static constexpr double PI_2 = static_cast(1.570796326794896619231321691639751442098584699687552910487L); + static constexpr double PI_4 = static_cast(0.785398163397448309615660845819875721049292349843776455243L); + namespace Matrix1234ComplexFloatBoost { + + // extern const BIG_COMPLEX_FLOAT TOLERANCE_LEVEL = BIG_COMPLEX_FLOAT(1e5, 1e5); // Initialization routine extern void Matrix1234Initializer(); @@ -59,6 +66,9 @@ namespace CFL_OBDD { extern CFLOBDD_COMPLEX_BIG MkPauliYMatrixInterleaved(unsigned int i); extern CFLOBDD_COMPLEX_BIG MkPauliZMatrixInterleaved(unsigned int i); extern CFLOBDD_COMPLEX_BIG MkSGateInterleaved(unsigned int i); + extern CFLOBDD_COMPLEX_BIG MkSdgGateInterleaved(unsigned int i); + extern CFLOBDD_COMPLEX_BIG MkTGateInterleaved(unsigned int i); + extern CFLOBDD_COMPLEX_BIG MkTdgGateInterleaved(unsigned int i); extern CFLOBDD_COMPLEX_BIG MkPhaseShiftGateInterleaved(unsigned int i, double theta); extern CFLOBDD_COMPLEX_BIG MkRestrictMatrix(unsigned int level, std::string s); @@ -77,6 +87,7 @@ namespace CFL_OBDD { // extern CFLOBDD_COMPLEX_BIG MatrixMultiply(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2); // Matrix multiplication extern CFLOBDD_COMPLEX_BIG MatrixMultiplyV4(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2); // Matrix multiplication extern CFLOBDD_COMPLEX_BIG MatrixMultiplyV4WithInfo(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2); // Matrix multiplication + extern CFLOBDD_COMPLEX_BIG MatrixTranspose(CFLOBDD_COMPLEX_BIG c); // Discrete Fourier Transform, and subroutines extern CFLOBDD_COMPLEX_BIG MkFourierMatrixInterleaved(unsigned int i); // Create representation of the DFT matrix @@ -103,6 +114,7 @@ namespace CFL_OBDD { extern CFLOBDD_COMPLEX_BIG MatrixShiftToAConnection(CFLOBDD_COMPLEX_BIG c); extern CFLOBDD_COMPLEX_BIG MatrixShiftToBConnection(CFLOBDD_COMPLEX_BIG c); extern CFLOBDD_COMPLEX_BIG KroneckerProduct2Vocs(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2); + extern CFLOBDD_COMPLEX_BIG ConjugateTranspose(CFLOBDD_COMPLEX_BIG c); } } diff --git a/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp b/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp index aabf771..0c782fe 100644 --- a/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp +++ b/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp @@ -30,6 +30,7 @@ #include #include #include +#include #include #include #include @@ -280,6 +281,77 @@ namespace CFL_OBDD { return v; } + CFLOBDDTopNodeComplexFloatBoostRefPtr MkSdgGateInterleavedTop(unsigned int i) + { + CFLOBDDTopNodeComplexFloatBoostRefPtr v; + CFLOBDDNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m; + + assert(i == 1); + + tempHandle = MkSGateInterleavedNode(i); + if (i >= 1){ + m.AddToEnd(1); + m.AddToEnd(0); + m.AddToEnd(BIG_COMPLEX_FLOAT(0, -1)); + } + // if (i >= 2) + // { + // m.AddToEnd(-1); + // } + // if (i >= 3) + // { + // m.AddToEnd(BIG_COMPLEX_FLOAT(0, -1)); + // } + m.Canonicalize(); + v = new CFLOBDDTopNodeComplexFloatBoost(tempHandle, m); + return v; + } + + CFLOBDDTopNodeComplexFloatBoostRefPtr MkTGateInterleavedTop(unsigned int i) + { + CFLOBDDTopNodeComplexFloatBoostRefPtr v; + CFLOBDDNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m; + + assert(i == 1); + + tempHandle = MkSGateInterleavedNode(i); + if (i >= 1){ + m.AddToEnd(1); + m.AddToEnd(0); + double cos_v = SQRT2_2; // boost::math::cos_pi(0.25); + double sin_v = SQRT2_2; // boost::math::sin_pi(0.25); + BIG_COMPLEX_FLOAT val(cos_v, sin_v); + m.AddToEnd(val); + } + m.Canonicalize(); + v = new CFLOBDDTopNodeComplexFloatBoost(tempHandle, m); + return v; + } + + CFLOBDDTopNodeComplexFloatBoostRefPtr MkTdgGateInterleavedTop(unsigned int i) + { + CFLOBDDTopNodeComplexFloatBoostRefPtr v; + CFLOBDDNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m; + + assert(i == 1); + + tempHandle = MkSGateInterleavedNode(i); + if (i >= 1){ + m.AddToEnd(1); + m.AddToEnd(0); + double cos_v = SQRT2_2; // boost::math::cos_pi(-0.25); + double sin_v = -SQRT2_2; // boost::math::sin_pi(-0.25); + BIG_COMPLEX_FLOAT val(cos_v, sin_v); + m.AddToEnd(val); + } + m.Canonicalize(); + v = new CFLOBDDTopNodeComplexFloatBoost(tempHandle, m); + return v; + } + CFLOBDDTopNodeComplexFloatBoostRefPtr MkSXGateTop(unsigned int i) { CFLOBDDTopNodeComplexFloatBoostRefPtr v; @@ -829,6 +901,51 @@ namespace CFL_OBDD { return new CFLOBDDTopNodeComplexFloatBoost(c, m); } + CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n) + { + std::unordered_map, + CFLOBDDNodeHandle::CFLOBDDNodeHandle_Hash> hashMap; + auto pc = MatrixTransposeNode(hashMap, *(n->rootConnection.entryPointHandle)); + CFLOBDDNodeHandle temp = pc.first; + ReductionMapHandle reductionMapHandle; + ComplexFloatBoostReturnMapHandle v; + for (unsigned int i = 0; i < pc.second.Size(); i++){ + reductionMapHandle.AddToEnd(i); + v.AddToEnd(n->rootConnection.returnMapHandle[pc.second[i]]); + } + reductionMapHandle.Canonicalize(); + v.Canonicalize(); + temp = temp.Reduce(reductionMapHandle, v.Size(), true); + return new CFLOBDDTopNodeComplexFloatBoost(temp, v); + return n; + } + + CFLOBDDTopNodeComplexFloatBoostRefPtr ConjugateTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c) + { + // std::unordered_map, + // CFLOBDDNodeHandle::CFLOBDDNodeHandle_Hash> hashMap; + // auto pc = MatrixTransposeNode(hashMap, *(c->rootConnection.entryPointHandle)); + // CFLOBDDNodeHandle temp = pc.first; + // ReductionMapHandle reductionMapHandle; + // ComplexFloatBoostReturnMapHandle v; + // for (unsigned int i = 0; i < pc.second.Size(); i++){ + // reductionMapHandle.AddToEnd(i); + // v.AddToEnd(c->rootConnection.returnMapHandle[pc.second[i]]); + // } + // reductionMapHandle.Canonicalize(); + // v.Canonicalize(); + // temp = temp.Reduce(reductionMapHandle, v.Size(), true); + // return new CFLOBDDTopNodeComplexFloatBoost(temp, v); + // return c; + ComplexFloatBoostReturnMapHandle m; + for (int i = 0; i < c->rootConnection.returnMapHandle.Size(); i++){ + auto val = c->rootConnection.returnMapHandle[i]; + m.AddToEnd(BIG_COMPLEX_FLOAT(val.real(), -val.imag())); + } + m.Canonicalize(); + return new CFLOBDDTopNodeComplexFloatBoost(*(c->rootConnection.entryPointHandle), m); + } + CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixMultiplyV4TopNode(CFLOBDDTopNodeComplexFloatBoostRefPtr c1, CFLOBDDTopNodeComplexFloatBoostRefPtr c2) { std::unordered_map hashMap; @@ -890,8 +1007,9 @@ namespace CFL_OBDD { CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixMultiplyV4WithInfoTopNode(CFLOBDDTopNodeComplexFloatBoostRefPtr c1, CFLOBDDTopNodeComplexFloatBoostRefPtr c2) { - if (c1->level >= 5) - clearMultMap(); + // if (c1->level >= 5) + // clearMultMap(); + const BIG_COMPLEX_FLOAT TOLERANCE_LEVEL = BIG_COMPLEX_FLOAT(1e10, 1e10); std::unordered_map hashMap; int c1_zero_index = -1, c2_zero_index = -1; c1_zero_index = c1->rootConnection.returnMapHandle.LookupInv(0); @@ -902,6 +1020,7 @@ namespace CFL_OBDD { ComplexFloatBoostReturnMapHandle v; boost::unordered_map reductionMap; ReductionMapHandle reductionMapHandle; + std::cout << std::setprecision(15); for (unsigned int i = 0; i < c->rootConnection.returnMapHandle.Size(); i++){ MatMultMapHandle r = c->rootConnection.returnMapHandle[i]; BIG_COMPLEX_FLOAT val = 0; @@ -910,31 +1029,29 @@ namespace CFL_OBDD { unsigned int index2 = j.first.second; if (index1 != -1 && index2 != -1){ auto factor = BIG_COMPLEX_FLOAT(j.second, 0); - // auto factor = j.second.convert_to(); - val = val + (factor * (c1->rootConnection.returnMapHandle[index1] * c2->rootConnection.returnMapHandle[index2])); + BIG_COMPLEX_FLOAT previous_val = val; + BIG_COMPLEX_FLOAT tmp_val = mul(c1->rootConnection.returnMapHandle[index1], c2->rootConnection.returnMapHandle[index2]); + BIG_COMPLEX_FLOAT added_val = mul(factor, tmp_val); + val = add(previous_val, added_val); } } + val = roundNearBy(val); + if (reductionMap.find(val) == reductionMap.end()){ v.AddToEnd(val); reductionMap.insert(std::make_pair(val, v.Size() - 1)); reductionMapHandle.AddToEnd(v.Size() - 1); + // std::cout << "New Value: " << val << " at " << v.Size() - 1 << std::endl; } else{ reductionMapHandle.AddToEnd(reductionMap[val]); + // std::cout << "Existing Value: " << val << " at " << reductionMap[val] << std::endl; } } - v.Canonicalize(); reductionMapHandle.Canonicalize(); CFLOBDDNodeHandle tempHandle = *(c->rootConnection.entryPointHandle); - // Perform reduction on tempHandle, with respect to the common elements that rmh maps together - /*ReductionMapHandle inducedReductionMapHandle; - FloatBoostReturnMapHandle inducedReturnMap; - v.InducedReductionAndReturnMap(inducedReductionMapHandle, inducedReturnMap); - CFLOBDDNodeHandle reduced_tempHandle = tempHandle.Reduce(inducedReductionMapHandle, inducedReturnMap.Size());*/ CFLOBDDNodeHandle reduced_tempHandle = tempHandle.Reduce(reductionMapHandle, v.Size(), true); - // Create and return CFLOBDDTopNode - //return(new CFLOBDDTopNodeFloatBoost(reduced_tempHandle, inducedReturnMap)); return(new CFLOBDDTopNodeComplexFloatBoost(reduced_tempHandle, v)); } } diff --git a/CFLOBDD/matrix1234_complex_float_boost_top_node.h b/CFLOBDD/matrix1234_complex_float_boost_top_node.h index c35e81d..19418f0 100644 --- a/CFLOBDD/matrix1234_complex_float_boost_top_node.h +++ b/CFLOBDD/matrix1234_complex_float_boost_top_node.h @@ -30,6 +30,7 @@ #include #include #include +#include #include "matrix1234_complex_float_boost.h" #include "return_map_T.h" namespace CFL_OBDD { @@ -63,6 +64,9 @@ namespace CFL_OBDD { extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkPauliYMatrixInterleavedTop(unsigned int i); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkPauliZMatrixInterleavedTop(unsigned int i); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSGateInterleavedTop(unsigned int i); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSdgGateInterleavedTop(unsigned int i); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkTGateInterleavedTop(unsigned int i); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkTdgGateInterleavedTop(unsigned int i); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkPhaseShiftGateInterleavedTop(unsigned int i, double theta); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkExchangeInterleavedTop(unsigned int i); // Representation of exchange matrix @@ -98,9 +102,106 @@ namespace CFL_OBDD { extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSXGateTop(unsigned int i); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSYGateTop(unsigned int i); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr ConjugateTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixShiftToAConnectionTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixShiftToBConnectionTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); + struct CustomComplexHash { + std::size_t operator()(const BIG_COMPLEX_FLOAT& p) const noexcept { + std::size_t seed = 0; + boost::hash_combine(seed, p.real()); + boost::hash_combine(seed, p.imag()); + return seed; + } + }; + + struct CustomComplexEqual { + bool operator()(const BIG_COMPLEX_FLOAT& lhs, const BIG_COMPLEX_FLOAT& rhs) const noexcept { + if (abs(lhs.real() - rhs.real()) < 1e-10 && + abs(lhs.imag() - rhs.imag()) < 1e-10) + return true; + return false; + } + }; + + inline BIG_COMPLEX_FLOAT isApproximatelyEqualTo(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b, double epsilon = 1e-10) { + return abs(a - b) < epsilon; + } + + template + T round_to(const T& val, int digits) { + using boost::multiprecision::round; // use multiprecision round + T factor = pow(T(10), digits); + return round(val * factor) / factor; + } + + inline BIG_COMPLEX_FLOAT roundNearBy(BIG_COMPLEX_FLOAT c) { + if (isApproximatelyEqualTo(c.real(), 0)) { + c = BIG_COMPLEX_FLOAT(0, c.imag()); + } else if (isApproximatelyEqualTo(c, 1)) { + c = BIG_COMPLEX_FLOAT(1, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(SQRT2_2, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(-SQRT2_2, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), 0.5)) { + c = BIG_COMPLEX_FLOAT(0.5, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -0.5)) { + c = BIG_COMPLEX_FLOAT(-0.5, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -1)) { + c = BIG_COMPLEX_FLOAT(-1, c.imag()); + } else { + c = BIG_COMPLEX_FLOAT(round_to(c.real(), 10), c.imag()); + } + + if (isApproximatelyEqualTo(c.imag(), 0)) { + c = BIG_COMPLEX_FLOAT(c.real(), 0); + } else if (isApproximatelyEqualTo(c.imag(), 1)) { + c = BIG_COMPLEX_FLOAT(c.real(), 1); + } else if (isApproximatelyEqualTo(c.imag(), -1)) { + c = BIG_COMPLEX_FLOAT(c.real(), -1); + } else if (isApproximatelyEqualTo(c.imag(), SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(c.real(), SQRT2_2); + } else if (isApproximatelyEqualTo(c.imag(), -SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(c.real(), -SQRT2_2); + } else if (isApproximatelyEqualTo(c.imag(), 0.5)) { + c = BIG_COMPLEX_FLOAT(c.real(), 0.5); + } else if (isApproximatelyEqualTo(c.imag(), -0.5)) { + c = BIG_COMPLEX_FLOAT(c.real(), -0.5); + } else { + c = BIG_COMPLEX_FLOAT(c.real(), round_to(c.imag(), 10)); + } + return c; + } + + inline BIG_COMPLEX_FLOAT div(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + if (b == 0) { + throw std::runtime_error("Division by zero in WeightedMatrix1234ComplexFloatBoostMul::div"); + } + auto c = a / b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT mul(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a * b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT add(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a + b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT sub(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a - b; + c = roundNearBy(c); + return c; + } + } } diff --git a/CFLOBDD/matrix1234_node.cpp b/CFLOBDD/matrix1234_node.cpp index d57a1b4..df062f1 100644 --- a/CFLOBDD/matrix1234_node.cpp +++ b/CFLOBDD/matrix1234_node.cpp @@ -3370,7 +3370,9 @@ namespace CFL_OBDD { } if (isFork) { - n->AConnection = Connection(CFLOBDDNodeHandle::CFLOBDDForkNodeHandle, commonly_used_return_maps[2]);//m01 + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); + n->AConnection = Connection(CFLOBDDNodeHandle::CFLOBDDForkNodeHandle, m01); n->numBConnections = 2; CFLOBDDReturnMapHandle m0; if (nhNode->BConnection[0].returnMapHandle.Lookup(0) == nhNode->BConnection[1].returnMapHandle.Lookup(0)){ @@ -3431,7 +3433,8 @@ namespace CFL_OBDD { } else { - n->AConnection = Connection(CFLOBDDNodeHandle::CFLOBDDDontCareNodeHandle, commonly_used_return_maps[0]);//m0 + CFLOBDDReturnMapHandle m0; m0.AddToEnd(0); m0.Canonicalize(); + n->AConnection = Connection(CFLOBDDNodeHandle::CFLOBDDDontCareNodeHandle, m0); n->numBConnections = 1; n->BConnection = new Connection[1]; n->BConnection[0] = Connection(CFLOBDDNodeHandle::CFLOBDDForkNodeHandle, nhNode->AConnection.returnMapHandle); @@ -5436,6 +5439,7 @@ namespace CFL_OBDD { CFLOBDDInternalNode* c1_internal = (CFLOBDDInternalNode *)c1.handleContents; CFLOBDDInternalNode* c2_internal = (CFLOBDDInternalNode *)c2.handleContents; + // std::cout << "MatrixMultiplyV4Node: AConnection level = " << c1.handleContents->level << std::endl; CFLOBDDTopNodeMatMultMapRefPtr aa = MatrixMultiplyV4Node(hashMap, *(c1_internal->AConnection.entryPointHandle), *(c2_internal->AConnection.entryPointHandle)); CFLOBDDReturnMapHandle mI; @@ -5449,15 +5453,21 @@ namespace CFL_OBDD { g->numExits = 0; //std::unordered_map mapFromHandleToIndex; std::unordered_map mapFromHandleToIndex; + // std::cout << "numBConnections: " << g->numBConnections << std::endl; + // for each B connection of g, compute the matrix multiplication of the corresponding for (unsigned int i = 0; i < g->numBConnections; i++){ MatMultMapHandle matmult_returnmap = aa->rootConnection.returnMapHandle[i]; CFLOBDDTopNodeMatMultMapRefPtr ans; bool first = true; // Consider Multiplication of M1 and M2 + // std::cout << "B Connection " << i << " has " << matmult_returnmap.mapContents->map.size() << " terms." << std::endl; + // matmult_returnmap.print(std::cout); for (auto &v : matmult_returnmap.mapContents->map){ unsigned int M1_index = v.first.first; unsigned int M2_index = v.first.second; //VAL_TYPE factor = v.second; + // std::cout << "Multiplying B connections: " << M1_index << " and " << M2_index << std::endl; + // Multiply B connections M1 and M2 CFLOBDDTopNodeMatMultMapRefPtr bb_old = MatrixMultiplyV4Node(hashMap, *(c1_internal->BConnection[M1_index].entryPointHandle), *(c2_internal->BConnection[M2_index].entryPointHandle)); @@ -5756,7 +5766,8 @@ namespace CFL_OBDD { hashMap.emplace(c2_zero_val_node, zm); } } - + + // std::cout << "MatrixMultiplyV4WithInfoNode: AConnection level = " << c1.handleContents->level << std::endl; CFLOBDDTopNodeMatMultMapRefPtr aa = MatrixMultiplyV4WithInfoNode(hashMap, *(c1_internal->AConnection.entryPointHandle), *(c2_internal->AConnection.entryPointHandle), a_zero_index_c1, a_zero_index_c2); CFLOBDDReturnMapHandle mI; @@ -5769,6 +5780,8 @@ namespace CFL_OBDD { g->BConnection = new Connection[g->numBConnections]; g->numExits = 0; //std::unordered_map mapFromHandleToIndex; + // std::cout << "numBConnections: " << g->numBConnections << std::endl; + // for each B connection of g, compute the matrix multiplication of the corresponding std::unordered_map mapFromHandleToIndex; for (unsigned int i = 0; i < g->numBConnections; i++){ MatMultMapHandle matmult_returnmap = aa->rootConnection.returnMapHandle[i]; @@ -5783,10 +5796,15 @@ namespace CFL_OBDD { } else{ // Consider Multiplication of M1 and M2 + // std::cout << "B Connection " << i << " has " << matmult_returnmap.mapContents->map.size() << " terms." << std::endl; + // matmult_returnmap.print(std::cout); + // for each term in the matrix multiplication map, multiply the corresponding B connections for (auto &v : matmult_returnmap.mapContents->map){ unsigned int M1_index = v.first.first; unsigned int M2_index = v.first.second; //VAL_TYPE factor = v.second; + // std::cout << "Multiplying B connections: " << M1_index << " and " << M2_index << std::endl; + // Multiply B connections M1 and M2 CFLOBDDTopNodeMatMultMapRefPtr bb_old = MatrixMultiplyV4WithInfoNode(hashMap, *(c1_internal->BConnection[M1_index].entryPointHandle), *(c2_internal->BConnection[M2_index].entryPointHandle), b_zero_indices_c1[M1_index], b_zero_indices_c2[M2_index]); @@ -5825,6 +5843,8 @@ namespace CFL_OBDD { else ans = bb; } + // std::cout << "Intermediate ans: after running M1_index = " << M1_index << " and M2_index = " << M2_index << std::endl; + // ans->rootConnection.print(std::cout); } } CFLOBDDReturnMapHandle ans_return_map; @@ -5843,6 +5863,10 @@ namespace CFL_OBDD { } } ans_return_map.Canonicalize(); + // std::cout << "Setting g->BConnection[" << i << "]" << std::endl; + // ans->rootConnection.entryPointHandle->print(std::cout); + // ans_return_map.print(std::cout); + // ans->rootConnection.returnMapHandle.print(std::cout); g->BConnection[i] = Connection(*(ans->rootConnection.entryPointHandle), ans_return_map); } @@ -5856,9 +5880,8 @@ namespace CFL_OBDD { #endif CFLOBDDNodeHandle gHandle(g); - //gHandle = gHandle.Reduce(reductionMapHandle, g_return_map.Size(), true); + gHandle = gHandle.Reduce(reductionMapHandle, g_return_map.Size(), true); CFLOBDDTopNodeMatMultMapRefPtr return_ans = new CFLOBDDTopNodeMatMultMap(gHandle, g_return_map); - //hashMap[mmp] = return_ans; matmult_hashMap_info[mmp] = return_ans; return return_ans; } diff --git a/CFLOBDD/tests_cfl.cpp b/CFLOBDD/tests_cfl.cpp index 2c7d6c8..71741b4 100644 --- a/CFLOBDD/tests_cfl.cpp +++ b/CFLOBDD/tests_cfl.cpp @@ -30,6 +30,7 @@ #include "matmult_map.h" // #include "matrix1234_node.h" #include "matrix1234_int.h" +#include "vector_complex_float_boost.h" #include "wmatrix1234_fb_mul.h" #include "weighted_cross_product.h" #include "wvector_fb_mul.h" @@ -1848,6 +1849,180 @@ void CFLTests::testSynBenchmark7_CFLOBDD(int size) std::cout << "Duration: " << duration.count() << " Memory: " << (nodeCount + edgeCount) << std::endl; } +// CFLOBDD_COMPLEX_BIG CreateGateF2(std::string indices, CFLOBDD_COMPLEX_BIG(*f)(unsigned int, unsigned int)) +// { +// if (indices.find('0') == std::string::npos) +// { +// unsigned int level = log2(indices.length() * 2); +// return f(level, 0); +// } +// else if (indices.find('1') == std::string::npos) +// { +// unsigned int level = log2(indices.length() * 2); +// return Matrix1234ComplexFloatBoost::MkIdRelationInterleaved(level); +// } +// else +// { +// auto F1 = CreateGateF2(indices.substr(0, indices.length()/2), f); +// auto F2 = CreateGateF2(indices.substr(indices.length()/2), f); +// return Matrix1234ComplexFloatBoost::KroneckerProduct2Vocs(F1, F2); +// } +// } + +// CFLOBDD_COMPLEX_BIG CreateGateF(std::string indices, CFLOBDD_COMPLEX_BIG(*f)(unsigned int)) +// { +// if (indices.find('0') == std::string::npos) +// { +// unsigned int level = log2(indices.length() * 2); +// return f(level); +// } +// else if (indices.find('1') == std::string::npos) +// { +// unsigned int level = log2(indices.length() * 2); +// return Matrix1234ComplexFloatBoost::MkIdRelationInterleaved(level); +// } +// else +// { +// auto F1 = CreateGateF(indices.substr(0, indices.length()/2), f); +// auto F2 = CreateGateF(indices.substr(indices.length()/2), f); +// return Matrix1234ComplexFloatBoost::KroneckerProduct2Vocs(F1, F2); +// } +// } + +WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL CreateGateF2_WCFLOBDD(std::string indices, WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(*f)(unsigned int, unsigned int, int)) +{ + if (indices.find('0') == std::string::npos) + { + unsigned int level = log2(indices.length() * 2); + return f(level, 0, 1); + } + else if (indices.find('1') == std::string::npos) + { + unsigned int level = log2(indices.length() * 2); + return WeightedMatrix1234ComplexFloatBoostMul::MkIdRelationInterleaved(level); + } + else + { + auto F1 = CreateGateF2_WCFLOBDD(indices.substr(0, indices.length()/2), f); + auto F2 = CreateGateF2_WCFLOBDD(indices.substr(indices.length()/2), f); + return WeightedMatrix1234ComplexFloatBoostMul::KroneckerProduct2Vocs(F1, F2); + } +} + +WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL CreateGateF_WCFLOBDD(std::string indices, WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(*f)(unsigned int, int, unsigned int)) +{ + if (indices.find('0') == std::string::npos) + { + unsigned int level = log2(indices.length() * 2); + return f(level, 1, 0); + } + else if (indices.find('1') == std::string::npos) + { + unsigned int level = log2(indices.length() * 2); + return WeightedMatrix1234ComplexFloatBoostMul::MkIdRelationInterleaved(level); + } + else + { + auto F1 = CreateGateF_WCFLOBDD(indices.substr(0, indices.length()/2), f); + auto F2 = CreateGateF_WCFLOBDD(indices.substr(indices.length()/2), f); + return WeightedMatrix1234ComplexFloatBoostMul::KroneckerProduct2Vocs(F1, F2); + } +} + +void CFLTests::testTmp() { + std::cout << std::setprecision(15); + auto state = CreateGateF2_WCFLOBDD("0010", WeightedVectorComplexFloatBoostMul::MkBasisVector); + // H (2) + auto IIHI = CreateGateF_WCFLOBDD("0010", WeightedMatrix1234ComplexFloatBoostMul::MkWalshInterleaved); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHI); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHI, state); + + // S (2) + auto IISdgI = CreateGateF_WCFLOBDD("0010", WeightedMatrix1234ComplexFloatBoostMul::MkSdgGate); + auto IISI = CreateGateF_WCFLOBDD("0010", WeightedMatrix1234ComplexFloatBoostMul::MkSGate); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IISdgI); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IISI, state); + + // CNOT (2, 3) + auto CNOT_2_3 = WeightedMatrix1234ComplexFloatBoostMul::MkCNOT(3, 4, 2, 3); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, CNOT_2_3); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(CNOT_2_3, state); + + // T (2) H (3) + auto IITdgI = CreateGateF_WCFLOBDD("0010", WeightedMatrix1234ComplexFloatBoostMul::MkTdgGate); + auto IIIH = CreateGateF_WCFLOBDD("0001", WeightedMatrix1234ComplexFloatBoostMul::MkWalshInterleaved); + auto IITI = CreateGateF_WCFLOBDD("0010", WeightedMatrix1234ComplexFloatBoostMul::MkTGate); + // auto IITdgH = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IITdgI, IIIH); + // auto IITH = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IITI, IIIH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IITdgH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IITH, state); + + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IITdgI); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IITI, state); + state.print(std::cout); + + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIIH); + state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIIH, state); + + return ; + + // // H (2) S (3) + // auto IIISdg = CreateGateF_WCFLOBDD("0001", WeightedMatrix1234ComplexFloatBoostMul::MkSdgGate); + // auto IIIS = CreateGateF_WCFLOBDD("0001", WeightedMatrix1234ComplexFloatBoostMul::MkSGate); + // auto IIHSdg = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHI, IIISdg); + // auto IIHS = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHI, IIIS); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHSdg); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHS, state); + + // // H (2) S (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHSdg); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHS, state); + + // // H (2) H (3) + // auto IIHH = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHI, IIIH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHH, state); + + // // S (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIISdg); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIIS, state); + + // // Sdg (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIIS); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIISdg, state); + + // // H (2) H (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHH, state); + + // // H (2) Sdg (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHS); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHSdg, state); + + // // H (2) Sdg (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHS); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHSdg, state); + + // // Tdg (2) H (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IITH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IITdgH, state); + + // // CNOT (2, 3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, CNOT_2_3); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(CNOT_2_3, state); + + // // Sdg (2) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IISI); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IISdgI, state); + + // // H (2) H (3) + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(state, IIHH); + // state = WeightedMatrix1234ComplexFloatBoostMul::MatrixMultiplyV4(IIHH, state); + + // state.print(std::cout); + +} + void CFLTests::InitModules() { @@ -2022,6 +2197,8 @@ bool CFLTests::runTests(const char *arg, int size, int seed, int a){ CFLTests::testSynBenchmark6_CFLOBDD(size); } else if (curTest == "testSyn7_CFL") { CFLTests::testSynBenchmark7_CFLOBDD(size); + } else if (curTest == "testTmp") { + CFLTests::testTmp(); } else { std::cout << "Unrecognized test name: " << curTest << std::endl; diff --git a/CFLOBDD/tests_cfl.h b/CFLOBDD/tests_cfl.h index 6bf9728..369cc90 100644 --- a/CFLOBDD/tests_cfl.h +++ b/CFLOBDD/tests_cfl.h @@ -54,6 +54,8 @@ class CFLTests static void testMatMul(int size); static void testQFT(int size, int seed); + static void testTmp(); + static void testWeightedOps(unsigned int size); static void testGHZAlgo_W(int size); static void testBVAlgo_W(int size, int seed); diff --git a/CFLOBDD/vector_complex_float_boost.cpp b/CFLOBDD/vector_complex_float_boost.cpp index 56c0ece..512bf45 100644 --- a/CFLOBDD/vector_complex_float_boost.cpp +++ b/CFLOBDD/vector_complex_float_boost.cpp @@ -103,6 +103,13 @@ namespace CFL_OBDD { { return GetPathCountTop(n.root, p); } + + long double ComputeL1Norm(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2) + { + auto c = m1 + (-1 * m2); + c.CountPaths(); + return ComputeNormTop(c.root); + } } } diff --git a/CFLOBDD/vector_complex_float_boost.h b/CFLOBDD/vector_complex_float_boost.h index ad6db3c..9b621e4 100644 --- a/CFLOBDD/vector_complex_float_boost.h +++ b/CFLOBDD/vector_complex_float_boost.h @@ -39,6 +39,7 @@ namespace CFL_OBDD { extern void VectorPrintColumnMajorInterleaved(CFLOBDD_COMPLEX_BIG c, std::ostream & out); extern long double getNonZeroProbability(CFLOBDD_COMPLEX_BIG n); extern unsigned long long int GetPathCount(CFLOBDD_COMPLEX_BIG n, long double p); + extern long double ComputeL1Norm(CFLOBDD_COMPLEX_BIG m1, CFLOBDD_COMPLEX_BIG m2); } } diff --git a/CFLOBDD/vector_complex_float_boost_top_node.cpp b/CFLOBDD/vector_complex_float_boost_top_node.cpp index 359121c..11097d7 100644 --- a/CFLOBDD/vector_complex_float_boost_top_node.cpp +++ b/CFLOBDD/vector_complex_float_boost_top_node.cpp @@ -187,6 +187,20 @@ namespace CFL_OBDD { return 0; } + long double ComputeNormTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n) + { + long double norm = 0.0; + for (unsigned int i = 0; i < n->rootConnection.returnMapHandle.Size(); i++) + { + if (n->rootConnection.returnMapHandle.Lookup(i) != 0){ + long double v = n->rootConnection.returnMapHandle.Lookup(i).real().convert_to(); + long double logNumPaths = n->rootConnection.entryPointHandle->handleContents->numPathsToExit[i]; + norm += std::pow(2, logNumPaths) * v * v; + } + } + return norm; + } + //#ifdef PATH_COUNTING_ENABLED std::string SamplingTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n, bool VocTwo) { diff --git a/CFLOBDD/vector_complex_float_boost_top_node.h b/CFLOBDD/vector_complex_float_boost_top_node.h index b88c93e..6209462 100644 --- a/CFLOBDD/vector_complex_float_boost_top_node.h +++ b/CFLOBDD/vector_complex_float_boost_top_node.h @@ -40,6 +40,7 @@ namespace CFL_OBDD { extern CFLOBDDTopNodeComplexFloatBoostRefPtr ConvertToDoubleTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); extern long double getNonZeroProbabilityTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n); extern unsigned long long int GetPathCountTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n, long double p); + extern long double ComputeNormTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n); //ifdef PATH_COUNTING_ENABLED extern std::string SamplingTop(CFLOBDDTopNodeComplexFloatBoostRefPtr n, bool voctwo = false); extern std::string SamplingV2Top(CFLOBDDTopNodeComplexFloatBoostRefPtr n); diff --git a/CFLOBDD/weighted_cflobdd_node_t.cpp b/CFLOBDD/weighted_cflobdd_node_t.cpp index 37776ec..e5ce103 100644 --- a/CFLOBDD/weighted_cflobdd_node_t.cpp +++ b/CFLOBDD/weighted_cflobdd_node_t.cpp @@ -79,14 +79,14 @@ void WeightedCFLOBDDNodeHandleT::InitNoDistinctionTable() n->AConnection.entryPointHandle = &(NoDistinctionNode[i-1]); m1.AddToEnd(0); m1.Canonicalize(); - n->AConnection.returnMapHandle = m1;//m1 + n->AConnection.returnMapHandle = m1;//m1 n->numBConnections = 1; n->BConnection = new WConnection[1]; n->BConnection[0].entryPointHandle = &(NoDistinctionNode[i-1]); m2.AddToEnd(0); m2.Canonicalize(); - n->BConnection[0].returnMapHandle = m2; + n->BConnection[0].returnMapHandle = m2; n->numExits = 1; //#ifdef PATH_COUNTING_ENABLED n->InstallPathCounts(); @@ -2027,7 +2027,7 @@ void WeightedCFLOBDDInternalNode::InstallWeightsOfPathsAsAmpsToExits() for (unsigned int i = 0; i < AConnection.entryPointHandle->handleContents->numExits; i++) { for (unsigned int j = 0; j < BConnection[i].entryPointHandle->handleContents->numExits; j++) { unsigned int k = BConnection[i].returnMapHandle.Lookup(j); - long double numPathsValue = AConnection.entryPointHandle->handleContents->numWeightsOfPathsAsAmpsToExit[i] * BConnection[i].entryPointHandle->handleContents->numWeightsOfPathsAsAmpsToExit[j]; + long double numPathsValue = AConnection.entryPointHandle->handleContents->numWeightsOfPathsAsAmpsToExit[i] * BConnection[i].entryPointHandle->handleContents->numWeightsOfPathsAsAmpsToExit[j]; // if (storingNumPathsToExit.find(k) == storingNumPathsToExit.end()){ // std::vector logOfPaths; // logOfPaths.push_back(numPathsValue); diff --git a/CFLOBDD/weighted_matmult_map.cpp b/CFLOBDD/weighted_matmult_map.cpp index 2cb545f..80b5121 100644 --- a/CFLOBDD/weighted_matmult_map.cpp +++ b/CFLOBDD/weighted_matmult_map.cpp @@ -139,8 +139,9 @@ bool WeightedMatMultMapBody::operator==(const WeightedMatMultMapBody &o) c else { for (auto m1_it = map.begin(), m2_it = o.map.begin(); m1_it != map.end() && m2_it != o.map.end(); m1_it++, m2_it++){ if ((m1_it->first.first != m2_it->first.first) || (m1_it->first.second != m2_it->first.second) || - ((m1_it->second != m2_it->second))) + ((m1_it->second != m2_it->second))) { return false; + } } } return true; @@ -277,6 +278,8 @@ template <> void WeightedMatMultMapHandle::Add(const INT_PAIR& p, BIG_COMPLEX_FLOAT& v) { assert(mapContents->refCount <= 1); + if (p.first == -1 || p.second == -1) + return ; auto it = mapContents->map.find(p); if (it == mapContents->map.end()){ mapContents->map.emplace(p, v); diff --git a/CFLOBDD/wmatrix1234_complex_fb_mul.cpp b/CFLOBDD/wmatrix1234_complex_fb_mul.cpp index 9c78c4e..c83b28a 100644 --- a/CFLOBDD/wmatrix1234_complex_fb_mul.cpp +++ b/CFLOBDD/wmatrix1234_complex_fb_mul.cpp @@ -49,6 +49,21 @@ namespace CFL_OBDD { return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkSGateTop(i, cflobdd_kind, offset)); } + WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkSdgGate(unsigned int i, int cflobdd_kind, unsigned int offset) + { + return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkSdgGateTop(i, cflobdd_kind, offset)); + } + + WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkTGate(unsigned int i, int cflobdd_kind, unsigned int offset) + { + return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkTGateTop(i, cflobdd_kind, offset)); + } + + WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkTdgGate(unsigned int i, int cflobdd_kind, unsigned int offset) + { + return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkTdgGateTop(i, cflobdd_kind, offset)); + } + WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkRestrictMatrix(unsigned int level, std::string s, int cflobdd_kind) { return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkRestrictTop(level, s, cflobdd_kind)); @@ -189,6 +204,11 @@ namespace CFL_OBDD { { return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(MkSYGateTop(level)); } + + WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL ConjugateTranspose(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL c) + { + return WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL(ConjugateTransposeTop(c.root)); + } } } diff --git a/CFLOBDD/wmatrix1234_complex_fb_mul.h b/CFLOBDD/wmatrix1234_complex_fb_mul.h index d196630..a03978e 100644 --- a/CFLOBDD/wmatrix1234_complex_fb_mul.h +++ b/CFLOBDD/wmatrix1234_complex_fb_mul.h @@ -15,8 +15,10 @@ namespace CFL_OBDD { typedef boost::multiprecision::cpp_complex_100 BIG_COMPLEX_FLOAT; // typedef double BIG_COMPLEX_FLOAT; typedef WEIGHTED_CFLOBDD_T> WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL; + // cflobdd_kind = 1 -> CFLOBDD, 0 -> BDD namespace WeightedMatrix1234ComplexFloatBoostMul { + // Initialization routine extern void Matrix1234Initializer(); @@ -26,6 +28,9 @@ namespace CFL_OBDD { extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkPauliYGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkPauliZGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkSGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkSdgGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkTGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkTdgGate(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkCNOTInterleaved(unsigned int i); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkExchangeInterleaved(unsigned int i); // Representation of exchange matrix extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkCNOT(unsigned int level, unsigned int n, long int controller, long int controlled, int cflobdd_kind = 1, unsigned int offset = 0); // Representation of CNOT matrix with index1 as controller and index2 as controlled bits @@ -53,6 +58,8 @@ namespace CFL_OBDD { extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkPhaseShiftGate(unsigned int i, double theta, int cflobdd_kind = 1, unsigned int offset = 0); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkCZGate(unsigned int i, long int c1, long int c2, double theta, int cflobdd_kind = 1, unsigned int offset = 0); extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL MkCSwapGate(unsigned int i, long int c1, long int x1, long int x2, int cflobdd_kind = 1, unsigned int offset = 0); + + extern WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL ConjugateTranspose(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL c); } } diff --git a/CFLOBDD/wmatrix1234_complex_fb_mul_node.cpp b/CFLOBDD/wmatrix1234_complex_fb_mul_node.cpp index b380ae9..074e617 100644 --- a/CFLOBDD/wmatrix1234_complex_fb_mul_node.cpp +++ b/CFLOBDD/wmatrix1234_complex_fb_mul_node.cpp @@ -320,39 +320,77 @@ namespace CFL_OBDD { { return std::make_pair(true, 0); } - else if (r1.mapContents->contains_zero_val == true || r2.mapContents->contains_zero_val == true) + bool are_equal = (r1 == r2); + return std::make_pair(are_equal, 1); + // else if (r1.mapContents->contains_zero_val == true || r2.mapContents->contains_zero_val == true) + // { + // return std::make_pair(false, 1); + // } + // else + // { + // auto it1 = r1.mapContents->map.begin(); + // auto it2 = r2.mapContents->map.begin(); + // for (it1 = r1.mapContents->map.begin(), it2 = r2.mapContents->map.begin(); + // it1 != r1.mapContents->map.end() && it2 != r2.mapContents->map.end(); + // it1++, it2++) + // { + // if (it1->first == it2->first) + // { + // if (it1->second == 0 || it2->second == 0) + // continue; + // if (first) + // factor = it2->second / it1->second; + // else{ + // BIG_COMPLEX_FLOAT f = it2->second / it1->second; + // if (f != factor) + // return std::make_pair(false, 1); + // } + // first = false; + // } + // else + // { + // return std::make_pair(false, 1); + // } + // } + // if (it1 != r1.mapContents->map.end() || it2 != r2.mapContents->map.end()) + // return std::make_pair(false, 1); + // return std::make_pair(true, factor); + // } + } + + std::pair, BIG_COMPLEX_FLOAT> factorize_if_needed(WeightedMatMultMapHandle& r) + { + BIG_COMPLEX_FLOAT factor = 1; + BIG_COMPLEX_FLOAT one = 1; + bool first = true; + WeightedMatMultMapHandle updated_r; + if (r.mapContents->contains_zero_val == true) { - return std::make_pair(false, 1); + return std::make_pair(r, 0); } else { - auto it1 = r1.mapContents->map.begin(); - auto it2 = r2.mapContents->map.begin(); - for (it1 = r1.mapContents->map.begin(), it2 = r2.mapContents->map.begin(); - it1 != r1.mapContents->map.end() && it2 != r2.mapContents->map.end(); - it1++, it2++) + bool first = true; + auto it = r.mapContents->map.begin(); + for (; it != r.mapContents->map.end(); + it++) { - if (it1->first == it2->first) - { - if (it1->second == 0 || it2->second == 0) - continue; - if (first) - factor = it2->second / it1->second; - else{ - BIG_COMPLEX_FLOAT f = it2->second / it1->second; - if (f != factor) - return std::make_pair(false, 1); + if (first) { + if (it->second == one) { + return std::make_pair(r, 1); } + factor = it->second; first = false; + updated_r.Add(it->first, one); } - else - { - return std::make_pair(false, 1); + else{ + auto weight = div(it->second, factor); + updated_r.Add(it->first, weight); } + } - if (it1 != r1.mapContents->map.end() || it2 != r2.mapContents->map.end()) - return std::make_pair(false, 1); - return std::make_pair(true, factor); + updated_r.Canonicalize(); + return std::make_pair(updated_r, factor); } } @@ -361,8 +399,26 @@ namespace CFL_OBDD { CFLOBDDMatMultMapHandle b1_m, CFLOBDDMatMultMapHandle b2_m, BIG_COMPLEX_FLOAT b_f1, BIG_COMPLEX_FLOAT b_f2) { WeightedPairProductMapHandle MapHandle; + // std::cout << "Adding Matrices\n"; + // std::cout << "Matrix 1 with bf1:" << b_f1 << "\n"; + // b1.print(std::cout); + // std::cout << "b1_m:\n"; + // for (size_t i = 0; i < b1_m.Size(); i++){ + // std::cout << "Entry " << i << ":\n"; + // b1_m.Lookup(i).print(std::cout); + // } + // std::cout << "Matrix 2 with bf2:" << b_f2 << "\n"; + // b2.print(std::cout); + // std::cout << "b2_m:\n"; + // for (size_t i = 0; i < b2_m.Size(); i++){ + // std::cout << "Entry " << i << ":\n"; + // b2_m.Lookup(i).print(std::cout); + // } auto tmp = PairProduct2(b1, b2, b_f1, b_f2, MapHandle); + // std::cout << "Resultant Matrix before reduction:\n"; + // tmp.print(std::cout); + CFLOBDDMatMultMapHandle returnMapHandle; boost::unordered_map reductionMap; @@ -390,8 +446,13 @@ namespace CFL_OBDD { val = v1 * c1 + v2 * c2; } val.Canonicalize(); + auto val_factored = factorize_if_needed(val); + val = val_factored.first; + BIG_COMPLEX_FLOAT factor_v = val_factored.second; + // std::cout << "Printing val_factored" << factor_v << "\n"; + // val.print(std::cout); bool found = false; - BIG_COMPLEX_FLOAT factor_c = 1; + BIG_COMPLEX_FLOAT factor_c = factor_v; unsigned int index = returnMapHandle.Size(); for (unsigned int k = 0; k < returnMapHandle.Size(); k++) { @@ -399,7 +460,7 @@ namespace CFL_OBDD { if (val_c.first == true) { found = true; - factor_c = val_c.second; + factor_c = factor_v; index = k; break; } @@ -409,9 +470,11 @@ namespace CFL_OBDD { // reductionMap.insert(std::make_pair(val.getHashCheck(), returnMapHandle.Size() - 1)); reductionMapHandle.AddToEnd(returnMapHandle.Size() - 1); if (val.mapContents->contains_zero_val == false) - valList.AddToEnd(1.0); + valList.AddToEnd(factor_c); else valList.AddToEnd(0.0); + // std::cout << "Added new entry to return map at " << returnMapHandle.Size() - 1 << "\n"; + // valList.print(std::cout); } else{ // reductionMapHandle.AddToEnd(reductionMap[val.getHashCheck()]); @@ -420,6 +483,8 @@ namespace CFL_OBDD { valList.AddToEnd(factor_c); else valList.AddToEnd(0.0); + // std::cout << "Reused existing entry in return map at " << index << "\n"; + // valList.print(std::cout); } iterator++; } @@ -429,6 +494,8 @@ namespace CFL_OBDD { valList.Canonicalize(); auto reduced_n = tmp.Reduce(reductionMapHandle, returnMapHandle.Size(), valList, false); BIG_COMPLEX_FLOAT factor = reduced_n.second; + // std::cout << "Resultant Matrix after reduction: " << factor << "\n"; + // reduced_n.first.print(std::cout); return std::make_tuple(reduced_n.first, returnMapHandle, factor); } @@ -550,16 +617,16 @@ namespace CFL_OBDD { int M2_numB = c2_internal->numBConnections; WeightedCFLOBDDComplexFloatBoostLeafNode* M1_b0 = (WeightedCFLOBDDComplexFloatBoostLeafNode *)c1_internal->BConnection[0].entryPointHandle->handleContents; WeightedCFLOBDDComplexFloatBoostLeafNode* M1_b1 = (WeightedCFLOBDDComplexFloatBoostLeafNode *)c1_internal->BConnection[M1_numB-1].entryPointHandle->handleContents; - a0 = a0 * M1_b0->lweight; - b0 = b0 * M1_b0->rweight; - c0 = c0 * M1_b1->lweight; - d0 = d0 * M1_b1->rweight; + a0 = mul(a0, M1_b0->lweight); + b0 = mul(b0, M1_b0->rweight); + c0 = mul(c0, M1_b1->lweight); + d0 = mul(d0, M1_b1->rweight); WeightedCFLOBDDComplexFloatBoostLeafNode* M2_b0 = (WeightedCFLOBDDComplexFloatBoostLeafNode *)c2_internal->BConnection[0].entryPointHandle->handleContents; WeightedCFLOBDDComplexFloatBoostLeafNode* M2_b1 = (WeightedCFLOBDDComplexFloatBoostLeafNode *)c2_internal->BConnection[M2_numB-1].entryPointHandle->handleContents; - a1 = a1 * M2_b0->lweight; - b1 = b1 * M2_b0->rweight; - c1 = c1 * M2_b1->lweight; - d1 = d1 * M2_b1->rweight; + a1 = mul(a1, M2_b0->lweight); + b1 = mul(b1, M2_b0->rweight); + c1 = mul(c1, M2_b1->lweight); + d1 = mul(d1, M2_b1->rweight); int M1_b0_numE = M1_b0->numExits; int M1_b1_numE = M1_b1->numExits; @@ -567,8 +634,8 @@ namespace CFL_OBDD { int M2_b1_numE = M2_b1->numExits; // v1 - BIG_COMPLEX_FLOAT a0a1 = a0 * a1; - BIG_COMPLEX_FLOAT b0c1 = b0 * c1; + BIG_COMPLEX_FLOAT a0a1 = mul(a0, a1); + BIG_COMPLEX_FLOAT b0c1 = mul(b0, c1); if (a0a1 == 0 && b0c1 == 0){ v1.Add(std::make_pair(-1,-1), zero); v1.mapContents->contains_zero_val = true; @@ -578,10 +645,17 @@ namespace CFL_OBDD { if (a0a1 != 0) v1.Add(std::make_pair(c1_internal->BConnection[0].returnMapHandle[0], c2_internal->BConnection[0].returnMapHandle[0]), a0a1); v1.Canonicalize(); + BIG_COMPLEX_FLOAT v1_factor = 1.0; + auto v1_factored_result = factorize_if_needed(v1); + v1 = v1_factored_result.first; + v1_factor = v1_factored_result.second; + // std::cout << "v1: " << v1_factor << std::endl; + // v1.print(std::cout); + // std::cout << std::endl; // v2 - BIG_COMPLEX_FLOAT a0b1 = a0 * b1; - BIG_COMPLEX_FLOAT b0d1 = b0 * d1; + BIG_COMPLEX_FLOAT a0b1 = mul(a0, b1); + BIG_COMPLEX_FLOAT b0d1 = mul(b0, d1); if (a0b1 == 0 && b0d1 == 0){ v2.Add(std::make_pair(-1,-1), zero); v2.mapContents->contains_zero_val = true; @@ -591,10 +665,17 @@ namespace CFL_OBDD { if (a0b1 != 0) v2.Add(std::make_pair(c1_internal->BConnection[0].returnMapHandle[0], c2_internal->BConnection[0].returnMapHandle[M2_b0_numE-1]), a0b1); v2.Canonicalize(); + BIG_COMPLEX_FLOAT v2_factor = 1.0; + auto v2_factored_result = factorize_if_needed(v2); + v2 = v2_factored_result.first; + v2_factor = v2_factored_result.second; + // std::cout << "v2: " << v2_factor << std::endl; + // v2.print(std::cout); + // std::cout << std::endl; // v3 - BIG_COMPLEX_FLOAT c0a1 = c0 * a1; - BIG_COMPLEX_FLOAT d0c1 = d0 * c1; + BIG_COMPLEX_FLOAT c0a1 = mul(c0, a1); + BIG_COMPLEX_FLOAT d0c1 = mul(d0, c1); if (c0a1 == 0 && d0c1 == 0){ v3.Add(std::make_pair(-1,-1), zero); v3.mapContents->contains_zero_val = true; @@ -604,10 +685,17 @@ namespace CFL_OBDD { if (c0a1 != 0) v3.Add(std::make_pair(c1_internal->BConnection[M1_numB-1].returnMapHandle[0], c2_internal->BConnection[0].returnMapHandle[0]), c0a1); v3.Canonicalize(); + BIG_COMPLEX_FLOAT v3_factor = 1.0; + auto v3_factored_result = factorize_if_needed(v3); + v3 = v3_factored_result.first; + v3_factor = v3_factored_result.second; + // std::cout << "v3: " << v3_factor << std::endl; + // v3.print(std::cout); + // std::cout << std::endl; // v4 - BIG_COMPLEX_FLOAT c0b1 = c0 * b1; - BIG_COMPLEX_FLOAT d0d1 = d0 * d1; + BIG_COMPLEX_FLOAT c0b1 = mul(c0, b1); + BIG_COMPLEX_FLOAT d0d1 = mul(d0, d1); if (c0b1 == 0 && d0d1 == 0){ v4.Add(std::make_pair(-1,-1), zero); v4.mapContents->contains_zero_val = true; @@ -617,168 +705,34 @@ namespace CFL_OBDD { if (c0b1 != 0) v4.Add(std::make_pair(c1_internal->BConnection[M1_numB-1].returnMapHandle[0], c2_internal->BConnection[0].returnMapHandle[M2_b0_numE-1]), c0b1); v4.Canonicalize(); - - // if (v1 == v3 && v2 == v4){ - // CFLOBDDReturnMapHandle m0; - // m0.AddToEnd(0); m0.Canonicalize(); - // if (v1 == v2){ - // g->numBConnections = 1; - // g->BConnection = new Connection[g->numBConnections]; - // // if (v1.mapContents->contains_zero_val == true){ - // // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m0); - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m0); - // // } - // // else { - // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m0); - // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m0); - // // } - // reductionMapHandle.AddToEnd(0); - // valList.AddToEnd(v1.mapContents->contains_zero_val ? 0 : 1); - // g_return_map.AddToEnd(v1); - // } - // else{ - // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m0); - // g->numBConnections = 1; - // g->BConnection = new Connection[g->numBConnections]; - // CFLOBDDReturnMapHandle m01; - // m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - // // if (v1.mapContents->contains_zero_val == true) - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01, m01); - // // else if (v2.mapContents->contains_zero_val == true) - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10, m01); - // // else - // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); - // reductionMapHandle.AddToEnd(0); - // reductionMapHandle.AddToEnd(1); - // valList.AddToEnd(v1.mapContents->contains_zero_val ? 0 : 1); - // valList.AddToEnd(v2.mapContents->contains_zero_val ? 0 : 1); - // g_return_map.AddToEnd(v1); - // g_return_map.AddToEnd(v2); - // } - // } - // else - // { - // CFLOBDDReturnMapHandle m01; - // m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); - // int la = 1, ra = 1; - // g->numBConnections = 2; - // g->BConnection = new Connection[2]; - // if (v1 == v2) - // { - // CFLOBDDReturnMapHandle m0; - // m0.AddToEnd(0); m0.Canonicalize(); - // // if (v1.mapContents->contains_zero_val == true){ - // // la = 0; - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m0); - // // } - // // else{ - // // la = 1; - // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m0); - // // } - // g_return_map.AddToEnd(v1); - // reductionMapHandle.AddToEnd(0); - // valList.AddToEnd(v1.mapContents->contains_zero_val ? 0 : 1); - // } - // else{ - // la = 1; - // // if (v1.mapContents->contains_zero_val == true) - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01, m01); - // // else if (v2.mapContents->contains_zero_val == true) - // // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10, m01); - // // else - // g->BConnection[0] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); - // g_return_map.AddToEnd(v1); - // g_return_map.AddToEnd(v2); - // reductionMapHandle.AddToEnd(0); - // reductionMapHandle.AddToEnd(1); - // valList.AddToEnd(v1.mapContents->contains_zero_val ? 0 : 1); - // valList.AddToEnd(v2.mapContents->contains_zero_val ? 0 : 1); - // } - // if (v3 == v4) - // { - // int k = 0; - // for (k = 0; k < g_return_map.Size(); k++) - // { - // if (g_return_map[k] == v3) - // break; - // } - // CFLOBDDReturnMapHandle m; - // m.AddToEnd(k); m.Canonicalize(); - // // if (v3.mapContents->contains_zero_val == true) - // // { - // // ra = 0; - // // g->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m); - // // } - // // else { - // // ra = 1; - // g->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m); - // // } - // if (k >= g_return_map.Size()){ - // g_return_map.AddToEnd(v3); - // reductionMapHandle.AddToEnd(k); - // valList.AddToEnd(v3.mapContents->contains_zero_val ? 0 : 1); - // } - // } - // else - // { - // ra = 1; - // int k1 = g_return_map.Size(), k2 = -1, k = 0; - // for (k = 0; k < g_return_map.Size(); k++) - // { - // if (g_return_map[k] == v3){ - // k1 = k; - // } - // else if (g_return_map[k] == v4){ - // k2 = k; - // } - // } - // if (k2 == -1){ - // k2 = std::max(k1, k-1) + 1; - // } - // CFLOBDDReturnMapHandle m1; - // m1.AddToEnd(k1); m1.AddToEnd(k2); m1.Canonicalize(); - // // if (v3.mapContents->contains_zero_val == true) - // // g->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01, m1); - // // else if (v4.mapContents->contains_zero_val == true) - // // g->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10, m1); - // // else - // g->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m1); - // if (k1 >= k){ - // g_return_map.AddToEnd(v3); - // reductionMapHandle.AddToEnd(k1); - // valList.AddToEnd(v3.mapContents->contains_zero_val ? 0 : 1); - // } - // if (k2 >= k){ - // g_return_map.AddToEnd(v4); - // reductionMapHandle.AddToEnd(k2); - // valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); - // } - // } - // // if (la == 1 && ra == 1) - // // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); - // // else if (la == 1 && ra == 0) - // // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10, m01); - // // else if (la == 0 && ra == 1) - // // g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01, m01); - // // assert(!(la == 0 && ra == 0)); - // } + BIG_COMPLEX_FLOAT v4_factor = 1.0; + auto v4_factored_result = factorize_if_needed(v4); + v4 = v4_factored_result.first; + v4_factor = v4_factored_result.second; + // std::cout << "v4: " << v4_factor << std::endl; + // v4.print(std::cout); + // std::cout << std::endl; g_return_map.AddToEnd(v1); valList.AddToEnd(v1.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(0); auto v2_c = compare(v1, v2); + Connection B1; + BIG_COMPLEX_FLOAT factor_b1 = 1.0; if (v2_c.first == true) { CFLOBDDReturnMapHandle m0; m0.AddToEnd(0); m0.Canonicalize(); - if (valList[0] == 0) + if (valList[0] == 0) { B1 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m0); + factor_b1 = 0.0; + } else { - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, v2_c.second)); + factor_b1 = v1_factor; + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, div(v2_factor, v1_factor))); B1 = Connection(x, m0); } } @@ -788,11 +742,22 @@ namespace CFL_OBDD { valList.AddToEnd(v2.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(1); CFLOBDDReturnMapHandle m01; m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - B1 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); + if (v1_factor == 0) { + B1 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m01); + factor_b1 = v2_factor; + } else if (v2_factor == 0) { + B1 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m01); + factor_b1 = v1_factor; + } else { + factor_b1 = v1_factor; + B1 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v2_factor, v1_factor)), m01); + } } auto v3_c = compare(v1, v3); Connection B2; + BIG_COMPLEX_FLOAT factor_b2 = 1.0; + // std::cout << "v3_c: " << v3_c.first << ", " << v3_c.second << std::endl; if (v3_c.first == true) { @@ -802,19 +767,22 @@ namespace CFL_OBDD { { CFLOBDDReturnMapHandle m0; m0.AddToEnd(0); m0.Canonicalize(); B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m0); + factor_b2 = 0.0; } else if (v2 == v4) { CFLOBDDReturnMapHandle m01; m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m01); + factor_b2 = v4_factor; } else { CFLOBDDReturnMapHandle m02; m02.AddToEnd(0); m02.AddToEnd(g_return_map.Size()); m02.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m02); + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m02); g_return_map.AddToEnd(v4); valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(g_return_map.Size()-1); + factor_b2 = v4_factor; } } else if (v4.mapContents->contains_zero_val == true) @@ -822,26 +790,31 @@ namespace CFL_OBDD { if (v2 == v4) { CFLOBDDReturnMapHandle m01; m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(v3_c.second, 0)); + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0)); B2 = Connection(x, m01); + factor_b2 = v3_factor; } else { CFLOBDDReturnMapHandle m02; m02.AddToEnd(0); m02.AddToEnd(g_return_map.Size()); m02.Canonicalize(); - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(v3_c.second, 0)); + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0)); B2 = Connection(x, m02); g_return_map.AddToEnd(v4); valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(g_return_map.Size()-1); + factor_b2 = v3_factor; } } else { auto v4_c = compare(v1, v4); + // std::cout << "v4_c: " << v4_c.first << ", " << v4_c.second << std::endl; + if (v4_c.first == true) { CFLOBDDReturnMapHandle m0; m0.AddToEnd(0); m0.Canonicalize(); - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(v3_c.second, v4_c.second)); + factor_b2 = v3_factor; + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, div(v4_factor, v3_factor))); B2 = Connection(x, m0); } else @@ -849,15 +822,17 @@ namespace CFL_OBDD { if (v2 == v4) { CFLOBDDReturnMapHandle m01; m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), m01); + factor_b2 = v3_factor; } else { CFLOBDDReturnMapHandle m02; m02.AddToEnd(0); m02.AddToEnd(g_return_map.Size()); m02.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m02); + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), m02); g_return_map.AddToEnd(v4); valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(g_return_map.Size()-1); + factor_b2 = v3_factor; } } } @@ -869,12 +844,27 @@ namespace CFL_OBDD { if (v2 == v4) { CFLOBDDReturnMapHandle mk; mk.AddToEnd(g_return_map.Size()-1); mk.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, mk); + if (v2_factor == 0) { + B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], mk); + factor_b2 = 0.0; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), mk); + factor_b2 = v3_factor; + } } else if (v1 == v4) { CFLOBDDReturnMapHandle mk; mk.AddToEnd(1); mk.AddToEnd(0); mk.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, mk); + if (v1_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), mk); + factor_b2 = v3_factor; + } else if (v2_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), mk); + factor_b2 = v4_factor; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), mk); + factor_b2 = v3_factor; + } } else { @@ -882,20 +872,32 @@ namespace CFL_OBDD { if (v4_c.first == true) { CFLOBDDReturnMapHandle m2; m2.AddToEnd(g_return_map.Size()-1); m2.Canonicalize(); - if (valList[valList.Size()-1] == 0) + if (valList[valList.Size()-1] == 0) { B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m2); + factor_b2 = 0.0; + } else { - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, v4_c.second)); + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, div(v4_factor, v3_factor))); B2 = Connection(x, m2); - } + factor_b2 = v3_factor; + } } else { valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(g_return_map.Size()); CFLOBDDReturnMapHandle m23; m23.AddToEnd(g_return_map.Size()-1); m23.AddToEnd(g_return_map.Size()); m23.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m23); + if (v3_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m23); + factor_b2 = v4_factor; + } else if (v4_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m23); + factor_b2 = v3_factor; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), m23); + factor_b2 = v3_factor; + } g_return_map.AddToEnd(v4); } } @@ -908,12 +910,24 @@ namespace CFL_OBDD { if (v1 == v4) { CFLOBDDReturnMapHandle m20; m20.AddToEnd(g_return_map.Size()-1); m20.AddToEnd(0); m20.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m20); + if (v1_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m20); + factor_b2 = v3_factor; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, v4_factor / v3_factor), m20); + factor_b2 = v3_factor; + } } else if (v2 == v4) { CFLOBDDReturnMapHandle m21; m21.AddToEnd(g_return_map.Size()-1); m21.AddToEnd(g_return_map.Size()-2); m21.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m21); + if (v2_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m21); + factor_b2 = v3_factor; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(v4_factor, v3_factor)), m21); + factor_b2 = v3_factor; + } } else { @@ -921,13 +935,16 @@ namespace CFL_OBDD { if (v4_c.first == true) { CFLOBDDReturnMapHandle m2; m2.AddToEnd(g_return_map.Size()-1); m2.Canonicalize(); - if (valList[valList.Size()-1] == 0) + if (valList[valList.Size()-1] == 0) { B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[0], m2); + factor_b2 = 0.0; + } else { - auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, v4_c.second)); + auto x = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, v4_factor / v3_factor)); B2 = Connection(x, m2); - } + factor_b2 = v3_factor; + } } else { @@ -935,17 +952,43 @@ namespace CFL_OBDD { valList.AddToEnd(v4.mapContents->contains_zero_val ? 0 : 1); reductionMapHandle.AddToEnd(g_return_map.Size()-1); CFLOBDDReturnMapHandle m23; m23.AddToEnd(g_return_map.Size()-2); m23.AddToEnd(g_return_map.Size()-1); m23.Canonicalize(); - B2 = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m23); + if (v3_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m23); + factor_b2 = v4_factor; + } else if (v4_factor == 0) { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m23); + factor_b2 = v3_factor; + } else { + B2 = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, v4_factor / v3_factor), m23); + factor_b2 = v3_factor; + } } } } } + // std::cout << "Curr c1 level 1: " << std::endl; + // c1_internal->print(std::cout); + // std::cout << "Curr c2 level 1: " << std::endl; + // c2_internal->print(std::cout); + + // std::cout << "B1: "; + // B1.print(std::cout); + // std::cout << " with factor " << factor_b1 << std::endl; + // std::cout << "B2: "; + // B2.print(std::cout); + // std::cout << " with factor " << factor_b2 << std::endl; + if (B1 == B2) { CFLOBDDReturnMapHandle m0; m0.AddToEnd(0); m0.Canonicalize(); - g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDDontCareNodeHandle, m0); + top_factor = mul(top_factor, factor_b1); + if (factor_b1 == 0) { + g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[1], m0); + } else { + g->AConnection = Connection(new WeightedCFLOBDDComplexFloatBoostDontCareNode(1, div(factor_b2, factor_b1)), m0); + } g->numBConnections = 1; g->BConnection = new Connection[1]; g->BConnection[0] = B1; @@ -953,15 +996,24 @@ namespace CFL_OBDD { else { CFLOBDDReturnMapHandle m01; m01.AddToEnd(0); m01.AddToEnd(1); m01.Canonicalize(); - g->AConnection = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle, m01); + if (factor_b1 == 0) { + g->AConnection = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(0, 1), m01); + top_factor = mul(top_factor, factor_b2); + } else if (factor_b2 == 0) { + g->AConnection = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, 0), m01); + top_factor = mul(top_factor, factor_b1); + } else { + g->AConnection = Connection(new WeightedCFLOBDDComplexFloatBoostForkNode(1, div(factor_b2, factor_b1)), m01); + top_factor = mul(top_factor, factor_b1); + } g->numBConnections = 2; g->BConnection = new Connection[2]; g->BConnection[0] = B1; g->BConnection[1] = B2; } - } - else{ + + } else { WeightedCFLOBDDComplexFloatBoostInternalNode* c1_internal = (WeightedCFLOBDDComplexFloatBoostInternalNode *)c1.handleContents; WeightedCFLOBDDComplexFloatBoostInternalNode* c2_internal = (WeightedCFLOBDDComplexFloatBoostInternalNode *)c2.handleContents; @@ -1064,12 +1116,12 @@ namespace CFL_OBDD { if (first){ ans = std::get<0>(bb_old); ans_matmult_map = new_bb_return; - ans_factor = v.second * std::get<2>(bb_old); + ans_factor = mul(v.second, std::get<2>(bb_old)); first = false; } else{ // TODO: change this - auto t = Add(ans, std::get<0>(bb_old), ans_matmult_map, new_bb_return, ans_factor, v.second * std::get<2>(bb_old)); + auto t = Add(ans, std::get<0>(bb_old), ans_matmult_map, new_bb_return, ans_factor, mul(v.second, std::get<2>(bb_old))); ans = std::get<0>(t); ans_matmult_map = std::get<1>(t); ans_factor = std::get<2>(t); @@ -1087,7 +1139,7 @@ namespace CFL_OBDD { } else { - ans_factor = ans_factor/factor; + ans_factor = div(ans_factor, factor); } // CFLOBDDReturnMapHandle ans_return_map; @@ -1163,7 +1215,7 @@ namespace CFL_OBDD { valList.AddToEnd(ans_factor); } - + g_return_map.Canonicalize(); reductionMapHandle.Canonicalize(); @@ -1188,7 +1240,9 @@ namespace CFL_OBDD { mA.Canonicalize(); g->AConnection = Connection(tmp.first, mA); WeightedCFLOBDDComplexFloatBoostMulNodeHandle gHandle(g); - ret = std::make_tuple(gHandle, g_return_map, tmp.second * top_factor * factor); + BIG_COMPLEX_FLOAT top_factor_final = mul(top_factor, factor); + top_factor_final = mul(top_factor_final, tmp.second); + ret = std::make_tuple(gHandle, g_return_map, top_factor_final); matmult_hash.insert(std::make_pair(mmp, ret)); for (unsigned int k = 0; k < nodeHandles.size(); k++) delete nodeHandles[k]; @@ -1211,13 +1265,13 @@ namespace CFL_OBDD { if (c1.handleContents->level == 1) { auto tmp = gHandle.Reduce(reductionMapHandle, g_return_map.Size(), valList, false); - ret = std::make_tuple(tmp.first, g_return_map, tmp.second * top_factor); + ret = std::make_tuple(tmp.first, g_return_map, mul(tmp.second, top_factor)); // ret = std::make_tuple(gHandle, g_return_map, top_factor); } else { auto tmp = gHandle.Reduce(reductionMapHandle, g_return_map.Size(), valList, true); - ret = std::make_tuple(tmp.first, g_return_map, tmp.second * top_factor); + ret = std::make_tuple(tmp.first, g_return_map, mul(tmp.second, top_factor)); // ret = std::make_tuple(gHandle, g_return_map, top_factor); } matmult_hash.insert(std::make_pair(mmp, ret)); @@ -2136,6 +2190,144 @@ namespace CFL_OBDD { return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(n); } + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkSdgGateNodeHelper(unsigned int i) + { + assert(i >= 1); + WeightedCFLOBDDComplexFloatBoostInternalNode *n = new WeightedCFLOBDDComplexFloatBoostInternalNode(i); + if (i == 1) { // Base case + + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(1, BIG_COMPLEX_FLOAT(0, -1))); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b0 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b1 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01; + CFLOBDDReturnMapHandle m10; + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + n->BConnection[0] = Connection(b0, m01); + n->BConnection[1] = Connection(b1, m10); + } + else { + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = MkSdgGateNodeHelper(i - 1); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + n->BConnection[0] = Connection(temp, m01); + CFLOBDDReturnMapHandle m1; + m1.AddToEnd(1); + m1.Canonicalize(); + n->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[i-1], m1); + } + n->numExits = 2; + #ifdef PATH_COUNTING_ENABLED + n->InstallPathCounts(); + #endif + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(n); + } + + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTGateNodeHelper(unsigned int i) + { + assert(i >= 1); + WeightedCFLOBDDComplexFloatBoostInternalNode *n = new WeightedCFLOBDDComplexFloatBoostInternalNode(i); + if (i == 1) { // Base case + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = + WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(1, BIG_COMPLEX_FLOAT(SQRT2_2, SQRT2_2))); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b0 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b1 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01; + CFLOBDDReturnMapHandle m10; + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + n->BConnection[0] = Connection(b0, m01); + n->BConnection[1] = Connection(b1, m10); + } + else { + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = MkTGateNodeHelper(i - 1); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + n->BConnection[0] = Connection(temp, m01); + CFLOBDDReturnMapHandle m1; + m1.AddToEnd(1); + m1.Canonicalize(); + n->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[i-1], m1); + } + n->numExits = 2; + #ifdef PATH_COUNTING_ENABLED + n->InstallPathCounts(); + #endif + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(n); + } + + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTdgGateNodeHelper(unsigned int i) + { + assert(i >= 1); + WeightedCFLOBDDComplexFloatBoostInternalNode *n = new WeightedCFLOBDDComplexFloatBoostInternalNode(i); + if (i == 1) { // Base case + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = + WeightedCFLOBDDComplexFloatBoostMulNodeHandle(new WeightedCFLOBDDComplexFloatBoostForkNode(1, BIG_COMPLEX_FLOAT(SQRT2_2, -SQRT2_2))); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b0 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle10; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle b1 = WeightedCFLOBDDComplexFloatBoostMulNodeHandle::CFLOBDDForkNodeHandle01; + CFLOBDDReturnMapHandle m10; + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + n->BConnection[0] = Connection(b0, m01); + n->BConnection[1] = Connection(b1, m10); + } + else { + WeightedCFLOBDDComplexFloatBoostMulNodeHandle temp = MkTdgGateNodeHelper(i - 1); + CFLOBDDReturnMapHandle m01; + m01.AddToEnd(0); + m01.AddToEnd(1); + m01.Canonicalize(); + n->AConnection = Connection(temp, m01); + + n->numBConnections = 2; + n->BConnection = new Connection[n->numBConnections]; + n->BConnection[0] = Connection(temp, m01); + CFLOBDDReturnMapHandle m1; + m1.AddToEnd(1); + m1.Canonicalize(); + n->BConnection[1] = Connection(WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[i-1], m1); + } + n->numExits = 2; + #ifdef PATH_COUNTING_ENABLED + n->InstallPathCounts(); + #endif + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(n); + } + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkPhaseShiftGateNodeHelper(unsigned int i, BIG_COMPLEX_FLOAT theta_val) { assert(i >= 1); @@ -3759,6 +3951,42 @@ namespace CFL_OBDD { } } + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkSdgGateNode(unsigned int level, int cflobdd_kind, unsigned int offset) + { + if (cflobdd_kind == 0) + { + abort(); + } + else + { + return MkSdgGateNodeHelper(level); + } + } + + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTGateNode(unsigned int level, int cflobdd_kind, unsigned int offset) + { + if (cflobdd_kind == 0) + { + abort(); + } + else + { + return MkTGateNodeHelper(level); + } + } + + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTdgGateNode(unsigned int level, int cflobdd_kind, unsigned int offset) + { + if (cflobdd_kind == 0) + { + abort(); + } + else + { + return MkTdgGateNodeHelper(level); + } + } + WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkiSwapGateNode(unsigned int level, long int index1, long int index2, int case_num, int cflobdd_kind, unsigned int offset) { if (cflobdd_kind == 0) @@ -3833,6 +4061,57 @@ namespace CFL_OBDD { } } + WeightedCFLOBDDComplexFloatBoostMulNodeHandle ConjugateTransposeNode(WeightedCFLOBDDComplexFloatBoostMulNodeHandle c) + { + unsigned int level = c.handleContents->level; + if (c == WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode_Ann[level] + || c == WeightedCFLOBDDComplexFloatBoostMulNodeHandle::NoDistinctionNode[level] + || c == WeightedCFLOBDDComplexFloatBoostMulNodeHandle::IdentityNode[level]) + return c; + + if (level == 0) + { + WeightedCFLOBDDComplexFloatBoostLeafNode* c1_internal = (WeightedCFLOBDDComplexFloatBoostLeafNode *)c.handleContents; + auto lweight = c1_internal->lweight; + auto rweight = c1_internal->rweight; + if (lweight.imag() == 0 && rweight.imag() == 0) + { + return c; + } + + auto new_lweight = lweight; + if (lweight.imag() != 0) + new_lweight = BIG_COMPLEX_FLOAT(lweight.real(), -lweight.imag()); + auto new_rweight = rweight; + if (rweight.imag() != 0) + new_rweight = BIG_COMPLEX_FLOAT(rweight.real(), -rweight.imag()); + if (c1_internal->numExits == 2) + { + WeightedCFLOBDDComplexFloatBoostForkNode* g = new WeightedCFLOBDDComplexFloatBoostForkNode(new_lweight, new_rweight); + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(g); + } else { + WeightedCFLOBDDComplexFloatBoostDontCareNode* g = new WeightedCFLOBDDComplexFloatBoostDontCareNode(new_lweight, new_rweight); + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(g); + } + } + else + { + WeightedCFLOBDDComplexFloatBoostInternalNode* c1_internal = (WeightedCFLOBDDComplexFloatBoostInternalNode *)c.handleContents; + WeightedCFLOBDDComplexFloatBoostInternalNode* g = new WeightedCFLOBDDComplexFloatBoostInternalNode(level); + g->numExits = c1_internal->numExits; + g->numBConnections = c1_internal->numBConnections; + auto aa = ConjugateTransposeNode(*c1_internal->AConnection.entryPointHandle); + g->AConnection = Connection(aa, c1_internal->AConnection.returnMapHandle); + g->BConnection = new Connection[g->numBConnections]; + for (unsigned int i = 0; i < g->numBConnections; i++) + { + auto bnode = ConjugateTransposeNode(*c1_internal->BConnection[i].entryPointHandle); + g->BConnection[i] = Connection(bnode, c1_internal->BConnection[i].returnMapHandle); + } + return WeightedCFLOBDDComplexFloatBoostMulNodeHandle(g); + } + } + } diff --git a/CFLOBDD/wmatrix1234_complex_fb_mul_node.h b/CFLOBDD/wmatrix1234_complex_fb_mul_node.h index 410b09c..7038f5d 100644 --- a/CFLOBDD/wmatrix1234_complex_fb_mul_node.h +++ b/CFLOBDD/wmatrix1234_complex_fb_mul_node.h @@ -10,6 +10,11 @@ namespace CFL_OBDD { + static constexpr double SQRT2_2 = static_cast(0.707106781186547524400844362104849039284835937688474036588L); + static constexpr double PI = static_cast(3.141592653589793238462643383279502884197169399375105820974L); + static constexpr double PI_2 = static_cast(1.570796326794896619231321691639751442098584699687552910487L); + static constexpr double PI_4 = static_cast(0.785398163397448309615660845819875721049292349843776455243L); + namespace WeightedMatrix1234ComplexFloatBoostMul { typedef boost::multiprecision::cpp_complex_100 BIG_COMPLEX_FLOAT; @@ -63,6 +68,9 @@ namespace CFL_OBDD { extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkPauliYGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkPauliZGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkSGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkSdgGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkTdgGateNode(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkCNOTInterleavedNode(unsigned int i); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkExchangeInterleavedNode(unsigned int i); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkCCPNode(std::unordered_map& cp_hashMap, unsigned int level, unsigned int n, long int controller1, long int controller2, long int controlled, BIG_COMPLEX_FLOAT theta); @@ -79,6 +87,7 @@ namespace CFL_OBDD { extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle MkCSwapGate2Node(unsigned int level, long int controller, long int i, long int j, int case_num, int cflobdd_kind = 1, unsigned int offset = 0); extern std::pair MkRestrictNode(unsigned int level, std::string s, int cflobdd_kind = 1); + extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle ConjugateTransposeNode(WeightedCFLOBDDComplexFloatBoostMulNodeHandle c); extern WeightedCFLOBDDComplexFloatBoostMulNodeHandle KroneckerProduct2VocsNode(WeightedCFLOBDDComplexFloatBoostMulNodeHandle m1, WeightedCFLOBDDComplexFloatBoostMulNodeHandle m2, int zero_index_m1, int zero_index_m2, bool rename = true); extern std::tuple @@ -93,6 +102,83 @@ namespace CFL_OBDD { // extern CFLOBDDTopNodeMatMultMapRefPtr MatrixMultiplyV4WithInfoNode( // std::unordered_map& hashMap, // WeightedCFLOBDDComplexFloatBoostMulNodeHandle c1, WeightedCFLOBDDComplexFloatBoostMulNodeHandle c2, int c1_zero_index, int c2_zero_index); + + inline BIG_COMPLEX_FLOAT isApproximatelyEqualTo(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b, double epsilon = 1e-10) { + return abs(a - b) < epsilon; + } + + template + T round_to(const T& val, int digits) { + using boost::multiprecision::round; // use multiprecision round + T factor = pow(T(10), digits); + return round(val * factor) / factor; + } + + inline BIG_COMPLEX_FLOAT roundNearBy(BIG_COMPLEX_FLOAT c) { + if (isApproximatelyEqualTo(c.real(), 0)) { + c = BIG_COMPLEX_FLOAT(0, c.imag()); + } else if (isApproximatelyEqualTo(c, 1)) { + c = BIG_COMPLEX_FLOAT(1, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(SQRT2_2, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(-SQRT2_2, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), 0.5)) { + c = BIG_COMPLEX_FLOAT(0.5, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -0.5)) { + c = BIG_COMPLEX_FLOAT(-0.5, c.imag()); + } else if (isApproximatelyEqualTo(c.real(), -1)) { + c = BIG_COMPLEX_FLOAT(-1, c.imag()); + } else { + c = BIG_COMPLEX_FLOAT(round_to(c.real(), 10), c.imag()); + } + + if (isApproximatelyEqualTo(c.imag(), 0)) { + c = BIG_COMPLEX_FLOAT(c.real(), 0); + } else if (isApproximatelyEqualTo(c.imag(), 1)) { + c = BIG_COMPLEX_FLOAT(c.real(), 1); + } else if (isApproximatelyEqualTo(c.imag(), -1)) { + c = BIG_COMPLEX_FLOAT(c.real(), -1); + } else if (isApproximatelyEqualTo(c.imag(), SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(c.real(), SQRT2_2); + } else if (isApproximatelyEqualTo(c.imag(), -SQRT2_2)) { + c = BIG_COMPLEX_FLOAT(c.real(), -SQRT2_2); + } else if (isApproximatelyEqualTo(c.imag(), 0.5)) { + c = BIG_COMPLEX_FLOAT(c.real(), 0.5); + } else if (isApproximatelyEqualTo(c.imag(), -0.5)) { + c = BIG_COMPLEX_FLOAT(c.real(), -0.5); + } else { + c = BIG_COMPLEX_FLOAT(c.real(), round_to(c.imag(), 10)); + } + return c; + } + + inline BIG_COMPLEX_FLOAT div(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + if (b == 0) { + throw std::runtime_error("Division by zero in WeightedMatrix1234ComplexFloatBoostMul::div"); + } + auto c = a / b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT mul(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a * b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT add(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a + b; + c = roundNearBy(c); + return c; + } + + inline BIG_COMPLEX_FLOAT sub(BIG_COMPLEX_FLOAT a, BIG_COMPLEX_FLOAT b) { + auto c = a - b; + c = roundNearBy(c); + return c; + } } } diff --git a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp index 4197b87..acb6852 100644 --- a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp +++ b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp @@ -98,6 +98,48 @@ namespace CFL_OBDD { return v; } + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkTdgGateTop(unsigned int i, int cflobdd_kind, unsigned int offset) + { + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr v; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m10; + + tempHandle = MkTdgGateNode(i, cflobdd_kind, offset); + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + v = new WeightedCFLOBDDTopNodeComplexFloatBoost(tempHandle, m10); + return v; + } + + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkSdgGateTop(unsigned int i, int cflobdd_kind, unsigned int offset) + { + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr v; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m10; + + tempHandle = MkSdgGateNode(i, cflobdd_kind, offset); + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + v = new WeightedCFLOBDDTopNodeComplexFloatBoost(tempHandle, m10); + return v; + } + + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkTGateTop(unsigned int i, int cflobdd_kind, unsigned int offset) + { + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr v; + WeightedCFLOBDDComplexFloatBoostMulNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m10; + + tempHandle = MkTGateNode(i, cflobdd_kind, offset); + m10.AddToEnd(1); + m10.AddToEnd(0); + m10.Canonicalize(); + v = new WeightedCFLOBDDTopNodeComplexFloatBoost(tempHandle, m10); + return v; + } + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkCNOTInterleavedTop(unsigned int i) { WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr v; @@ -164,12 +206,12 @@ namespace CFL_OBDD { tempHandle = MkWalshInterleavedNode(i, cflobdd_kind, offset); BIG_COMPLEX_FLOAT val; if (cflobdd_kind == 0) - val = boost::multiprecision::pow(BIG_COMPLEX_FLOAT(sqrt(2)), i/2).convert_to(); + val = 1.0/boost::multiprecision::pow(BIG_COMPLEX_FLOAT(sqrt(2)), i/2).convert_to(); else - val = boost::multiprecision::pow(BIG_COMPLEX_FLOAT(sqrt(2)), i-1).convert_to(); + val = boost::multiprecision::pow(BIG_COMPLEX_FLOAT(SQRT2_2), boost::multiprecision::pow(BIG_COMPLEX_FLOAT(2), i-1)).convert_to(); m.AddToEnd(1); m.Canonicalize(); - v = new WeightedCFLOBDDTopNodeComplexFloatBoost(tempHandle, m, 1.0/val); + v = new WeightedCFLOBDDTopNodeComplexFloatBoost(tempHandle, m, val); return v; } @@ -363,8 +405,11 @@ namespace CFL_OBDD { long int index1 = j.first.first; long int index2 = j.first.second; BIG_COMPLEX_FLOAT factor(j.second); - if (!(index1 == -1 && index2 == -1)) - val = val + (factor * (c1->rootConnection.returnMapHandle[index1] * c2->rootConnection.returnMapHandle[index2])); + if (!(index1 == -1 && index2 == -1)) { + auto tmp_val = mul(c1->rootConnection.returnMapHandle[index1], c2->rootConnection.returnMapHandle[index2]); + tmp_val = mul(tmp_val, factor); + val = add(val, tmp_val); + } } if (first && val != 0) { @@ -374,7 +419,7 @@ namespace CFL_OBDD { } else { - val = val/common_f; + val = div(val, common_f); } BIG_COMPLEX_FLOAT val_to_check = (val == 0)? 0 : 1; if (reductionMap.find(val_to_check) == reductionMap.end()){ @@ -396,7 +441,11 @@ namespace CFL_OBDD { auto g = n.Reduce(reductionMapHandle, v.Size(), valList, true); // Create and return CFLOBDDTopNode //return(new CFLOBDDTopNodeFloatBoost(reduced_tempHandle, inducedReturnMap)); - return(new WeightedCFLOBDDTopNodeComplexFloatBoost(g.first, v, g.second * common_f * c1->rootConnection.factor * c2->rootConnection.factor * n_factor)); + auto top_factor = mul(common_f, n_factor); + top_factor = mul(top_factor, c1->rootConnection.factor); + top_factor = mul(top_factor, c2->rootConnection.factor); + top_factor = mul(top_factor, g.second); + return(new WeightedCFLOBDDTopNodeComplexFloatBoost(g.first, v, top_factor)); } WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkCNOTTopNode(unsigned int level, unsigned int n, long int controller, long int controlled, int cflobdd_kind, unsigned int offset) @@ -629,6 +678,16 @@ namespace CFL_OBDD { return v; } + WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr ConjugateTransposeTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr c) + { + // only does conjugate. doesn't do transpose + WeightedCFLOBDDComplexFloatBoostMulNodeHandle ct = ConjugateTransposeNode(*(c->rootConnection.entryPointHandle)); + BIG_COMPLEX_FLOAT conj_factor = c->rootConnection.factor; + if (conj_factor.imag() != 0) + conj_factor = BIG_COMPLEX_FLOAT(conj_factor.real(), -conj_factor.imag()); + return new WeightedCFLOBDDTopNodeComplexFloatBoost(ct, c->rootConnection.returnMapHandle, conj_factor); + } + } } diff --git a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.h b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.h index b0e388d..b85d589 100644 --- a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.h +++ b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.h @@ -36,6 +36,9 @@ namespace CFL_OBDD { extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkPauliYGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkPauliZGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkSGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkSdgGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkTGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); + extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkTdgGateTop(unsigned int i, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkCNOTInterleavedTop(unsigned int i); extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkExchangeInterleavedTop(unsigned int i); // Representation of exchange matrix extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkCNOTTopNode(unsigned int level, unsigned int n, long int controller, long int controlled, int cflobdd_kind = 1, unsigned int offset = 0); @@ -43,6 +46,7 @@ namespace CFL_OBDD { extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkCCNOTTop(unsigned int level, long int controller1, long int controller2, long int controlled, int cflobdd_kind = 1, unsigned int offset = 0); extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr MkRestrictTop(unsigned int level, std::string s, int cflobdd_kind = 1); + extern WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr ConjugateTransposeTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr c); extern void MatrixPrintRowMajorTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr n, std::ostream & out); extern void MatrixPrintRowMajorInterleavedTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr n, std::ostream & out); diff --git a/CFLOBDD/wvector_complex_fb_mul.cpp b/CFLOBDD/wvector_complex_fb_mul.cpp index a3f2f5f..c56405a 100644 --- a/CFLOBDD/wvector_complex_fb_mul.cpp +++ b/CFLOBDD/wvector_complex_fb_mul.cpp @@ -80,6 +80,14 @@ namespace CFL_OBDD { { return getNonZeroProbabilityTop(n.root); } + + long double ComputeL1Norm(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL m1, WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL m2) + { + BIG_COMPLEX_FLOAT minus_one (-1, 0); + auto c = m1 + (minus_one * m2); + c.ComputeWeightOfPathsAsAmpsToExits(); + return ComputeNormTop(c.root); + } } } diff --git a/CFLOBDD/wvector_complex_fb_mul.h b/CFLOBDD/wvector_complex_fb_mul.h index 981174d..1ec8228 100644 --- a/CFLOBDD/wvector_complex_fb_mul.h +++ b/CFLOBDD/wvector_complex_fb_mul.h @@ -32,6 +32,8 @@ namespace CFL_OBDD { extern void VectorPrintColumnMajorInterleaved(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL c, std::ostream & out); extern void PrintVector(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL c, std::ostream & out, unsigned int vars_count); extern long double getNonZeroProbability(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL n); + + extern long double ComputeL1Norm(WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL m1, WEIGHTED_CFLOBDD_COMPLEX_FLOAT_BOOST_MUL m2); } } diff --git a/CFLOBDD/wvector_complex_fb_mul_node.cpp b/CFLOBDD/wvector_complex_fb_mul_node.cpp index a9dbb5c..1839dd3 100644 --- a/CFLOBDD/wvector_complex_fb_mul_node.cpp +++ b/CFLOBDD/wvector_complex_fb_mul_node.cpp @@ -17,6 +17,7 @@ namespace CFL_OBDD { std::vector> commonly_used_return_maps_vector;// m0, m1, m01, m10 + void InitReturnMapVectorHandles(){ ReturnMapHandle m0, m1, m01, m10; m0.AddToEnd(0); diff --git a/CFLOBDD/wvector_top_node_complex_fb_mul.cpp b/CFLOBDD/wvector_top_node_complex_fb_mul.cpp index fe3b64c..3fb6fd5 100644 --- a/CFLOBDD/wvector_top_node_complex_fb_mul.cpp +++ b/CFLOBDD/wvector_top_node_complex_fb_mul.cpp @@ -231,6 +231,18 @@ namespace CFL_OBDD { long double prob = getNonZeroProbabilityNode(*(n->rootConnection.entryPointHandle)); return prob / factor.convert_to(); } + + long double ComputeNormTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr c) + { + if (c->rootConnection.factor == 0) + return 0; + int index = c->rootConnection.returnMapHandle.LookupInv(1); + auto norm = c->rootConnection.entryPointHandle->handleContents->numWeightsOfPathsAsAmpsToExit[index]; + long double real = c->rootConnection.factor.real().convert_to(); + long double imag = c->rootConnection.factor.imag().convert_to(); + long double factor = real * real + imag * imag; + return norm * factor; + } } } diff --git a/CFLOBDD/wvector_top_node_complex_fb_mul.h b/CFLOBDD/wvector_top_node_complex_fb_mul.h index 242455f..2f10112 100644 --- a/CFLOBDD/wvector_top_node_complex_fb_mul.h +++ b/CFLOBDD/wvector_top_node_complex_fb_mul.h @@ -38,6 +38,7 @@ namespace CFL_OBDD { extern void VectorPrintColumnMajorInterleavedTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr n, std::ostream & out); extern void PrintVectorTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr n, std::ostream & out, unsigned int vars_count); extern long double getNonZeroProbabilityTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr n); + extern long double ComputeNormTop(WeightedCFLOBDDTopNodeComplexFloatBoostRefPtr c); } } From fffd78d88fbd7be8cdfb77b9226295d1d4f98849 Mon Sep 17 00:00:00 2001 From: Meghana Sistla Date: Mon, 29 Dec 2025 23:27:14 +0000 Subject: [PATCH 2/3] Chnages for equivalence checking --- CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp index acb6852..dd9f86e 100644 --- a/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp +++ b/CFLOBDD/wmatrix1234_top_node_complex_fb_mul.cpp @@ -394,8 +394,21 @@ namespace CFL_OBDD { ReductionMapHandle reductionMapHandle; WeightedValuesListHandle valList; WeightedCFLOBDDComplexFloatBoostMulNodeHandle n = std::get<0>(c); + // std::cout << "Printing n:" << std::endl; + // n.print(std::cout); CFLOBDDMatMultMapHandle n_return = std::get<1>(c); + // std::cout << "Printing n_return:" << std::endl; + // for (unsigned int i = 0; i < n_return.Size(); i++){ + // WeightedMatMultMapHandle r = n_return[i]; + // std::cout << "n_return[" << i << "]:" << std::endl; + // r.print(std::cout); + // } BIG_COMPLEX_FLOAT n_factor = std::get<2>(c); + // std::cout << "Printing n_factor: " << n_factor << std::endl; + // std::cout << "Printing c1 return map:" << std::endl; + // c1->rootConnection.returnMapHandle.print(std::cout); + // std::cout << "Printing c2 return map:" << std::endl; + // c2->rootConnection.returnMapHandle.print(std::cout); bool first = true; BIG_COMPLEX_FLOAT common_f = 1.0; for (unsigned int i = 0; i < n_return.Size(); i++){ @@ -436,8 +449,9 @@ namespace CFL_OBDD { v.Canonicalize(); reductionMapHandle.Canonicalize(); valList.Canonicalize(); - // n.print(std::cout); // std::cout << valList << std::endl; + // std::cout << "Printing reduction map:" << std::endl; + // reductionMapHandle.print(std::cout); auto g = n.Reduce(reductionMapHandle, v.Size(), valList, true); // Create and return CFLOBDDTopNode //return(new CFLOBDDTopNodeFloatBoost(reduced_tempHandle, inducedReturnMap)); From 9120f55152fc096d9a8dc69a50aba3943490fe91 Mon Sep 17 00:00:00 2001 From: Meghana Sistla Date: Sun, 18 Jan 2026 17:59:05 +0000 Subject: [PATCH 3/3] Adding RZ Gate --- CFLOBDD/matrix1234_complex_float_boost.cpp | 6 ++++++ CFLOBDD/matrix1234_complex_float_boost.h | 1 + ...atrix1234_complex_float_boost_top_node.cpp | 21 +++++++++++++++++++ .../matrix1234_complex_float_boost_top_node.h | 1 + 4 files changed, 29 insertions(+) diff --git a/CFLOBDD/matrix1234_complex_float_boost.cpp b/CFLOBDD/matrix1234_complex_float_boost.cpp index 0b52379..2d075de 100644 --- a/CFLOBDD/matrix1234_complex_float_boost.cpp +++ b/CFLOBDD/matrix1234_complex_float_boost.cpp @@ -133,6 +133,12 @@ namespace CFL_OBDD { } + CFLOBDD_COMPLEX_BIG MkRZGate(unsigned int i, double theta) + { + return CFLOBDD_COMPLEX_BIG(MkRZGateTop(i, theta)); + + } + // Create representation of the Inverse Reed-Muller matrix IRM(2**(i-1)) // [i.e., a matrix of size 2**(2**(i-1))) x 2**(2**(i-1)))] // with interleaved indexing of components: that is, input diff --git a/CFLOBDD/matrix1234_complex_float_boost.h b/CFLOBDD/matrix1234_complex_float_boost.h index e37653d..57da35e 100644 --- a/CFLOBDD/matrix1234_complex_float_boost.h +++ b/CFLOBDD/matrix1234_complex_float_boost.h @@ -109,6 +109,7 @@ namespace CFL_OBDD { extern CFLOBDD_COMPLEX_BIG MkCCP(unsigned int level, unsigned int n, long int controller1, long int controller2, long int controlled, double theta); extern CFLOBDD_COMPLEX_BIG MkSXGate(unsigned int i); extern CFLOBDD_COMPLEX_BIG MkSYGate(unsigned int i); + extern CFLOBDD_COMPLEX_BIG MkRZGate(unsigned int i, double theta); extern CFLOBDD_COMPLEX_BIG MatrixShiftToAConnection(CFLOBDD_COMPLEX_BIG c); diff --git a/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp b/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp index 0c782fe..a7a75e2 100644 --- a/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp +++ b/CFLOBDD/matrix1234_complex_float_boost_top_node.cpp @@ -352,6 +352,27 @@ namespace CFL_OBDD { return v; } + CFLOBDDTopNodeComplexFloatBoostRefPtr MkRZGateTop(unsigned int i, double theta) + { + CFLOBDDTopNodeComplexFloatBoostRefPtr v; + CFLOBDDNodeHandle tempHandle; + ComplexFloatBoostReturnMapHandle m; + + assert(i == 1); + + double cos_v = boost::math::cos_pi(theta); + double sin_v = boost::math::sin_pi(theta); + BIG_COMPLEX_FLOAT val1(cos_v, sin_v); + BIG_COMPLEX_FLOAT val2(-cos_v, -sin_v); + tempHandle = MkSGateInterleavedNode(i); + m.AddToEnd(val1); + m.AddToEnd(0); + m.AddToEnd(val2); + m.Canonicalize(); + v = new CFLOBDDTopNodeComplexFloatBoost(tempHandle, m); + return v; + } + CFLOBDDTopNodeComplexFloatBoostRefPtr MkSXGateTop(unsigned int i) { CFLOBDDTopNodeComplexFloatBoostRefPtr v; diff --git a/CFLOBDD/matrix1234_complex_float_boost_top_node.h b/CFLOBDD/matrix1234_complex_float_boost_top_node.h index 19418f0..cfd873d 100644 --- a/CFLOBDD/matrix1234_complex_float_boost_top_node.h +++ b/CFLOBDD/matrix1234_complex_float_boost_top_node.h @@ -101,6 +101,7 @@ namespace CFL_OBDD { extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkCCPTopNode(unsigned int level, unsigned int n, long int controller1, long int controller2, long int controlled, double theta); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSXGateTop(unsigned int i); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkSYGateTop(unsigned int i); + extern CFLOBDDTopNodeComplexFloatBoostRefPtr MkRZGateTop(unsigned int i, double theta); extern CFLOBDDTopNodeComplexFloatBoostRefPtr ConjugateTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c); extern CFLOBDDTopNodeComplexFloatBoostRefPtr MatrixTransposeTop(CFLOBDDTopNodeComplexFloatBoostRefPtr c);