Add parallel option

.github/workflows/ci.yml

@@ -2,38 +2,46 @@ name: QuantumSim CI
 
 on:
   push:
-    branches: [master]
+    branches: [main]
   pull_request:
-    branches: [master]
+    branches: [main]
 
 jobs:
   build-ubuntu:
-    runs-on: ubuntu-latest
+    runs-on: ${{ matrix.os }}
+
+    strategy:
+      matrix:
+        os: [ubuntu-latest, ubuntu-24.04, ubuntu-22.04]
 
     steps:
       - uses: actions/checkout@v2
+      - name: Install OpenMP
+        run: |
+          sudo apt-get update
+          sudo apt-get install -y libomp-dev
       - name: make check
         run: make check
       - name: make
         run: make
 
   build-macos:
-    runs-on: macos-latest
-
-    steps:
-      - uses: actions/checkout@v2
-      - name: install libomp
-        run: brew install libomp
-      - name: make check
-        run: make check
-      - name: make
-        run: make
+    runs-on: ${{ matrix.os }}
 
-  build-macos-no-openmp:
-    runs-on: macos-latest
+    strategy:
+      matrix:
+        os: [macos-latest, macos-15, macos-14, macos-13, macos-12, macos-11]
 
     steps:
       - uses: actions/checkout@v2
+      - name: Install OpenMP
+        run: brew install libomp
+      - name: Set environment variables for libomp
+        run: |
+          export LDFLAGS="-L$(brew --prefix libomp)/lib"
+          export CPPFLAGS="-I$(brew --prefix libomp)/include"
+          echo "LDFLAGS=$LDFLAGS" >> $GITHUB_ENV
+          echo "CPPFLAGS=$CPPFLAGS" >> $GITHUB_ENV
       - name: make check
         run: make check
       - name: make

.gitignore

@@ -1,24 +1,11 @@
+.DS_STORE
+
 # vscode directory
 .vscode/
 
 # python cache
 __pycache__/
 
-.DS_STORE
-
-# data
-data/
-
-# plotting script
-experiments.py
-
-# images, figures etc.
-figs/
-*.eps
-
-# transfer script
-transfer
-
 # object files
 *.o
 

Makefile

@@ -7,7 +7,7 @@ CXX = g++
 # set linker flag for OpenMP based on OS
 ifeq ($(OS_NAME), darwin)
 OPENMP_LINKER_FLAG = -lomp
-OPENMP_FLAGS = -Xpreprocessor -fopenmp
+OPENMP_FLAGS = -Xclang -fopenmp
 endif
 ifeq ($(OS_NAME), linux)
 OPENMP_LINKER_FLAG = -lgomp
@@ -94,12 +94,12 @@ $(TARGET): $(TARGET).o $(QUBITLAYER).o $(EXAMPLES).o
 
 $(TARGET).o: $(TARGET).cpp $(TARGET_DEPS) $(QLAYER_DEPS)
 	@printf "%b" "$(BLUE)$(COM_STRING) $(NO_COLOR)$(@)                             				"
-	@$(CXX) $(CXXFLAGS) -c $(TARGET).cpp -o $(TARGET).o
+	@$(CXX) $(CXXFLAGS) $(OPENMP_FLAGS) -c $(TARGET).cpp -o $(TARGET).o
 	@printf "%b" "$(GREEN)$(OK_STRING)$(NO_COLOR)\n"
 
 $(QUBITLAYER).o: $(QUBITLAYER).cpp $(TARGET_DEPS) $(QLAYER_DEPS)
 	@printf "%b" "$(BLUE)$(COM_STRING) $(NO_COLOR)$(@)                       				"
-	@$(CXX) $(CXXFLAGS) -c $(QUBITLAYER).cpp -o $(QUBITLAYER).o
+	@$(CXX) $(CXXFLAGS) $(OPENMP_FLAGS) -c $(QUBITLAYER).cpp -o $(QUBITLAYER).o
 	@printf "%b" "$(GREEN)$(OK_STRING)$(NO_COLOR)\n"
 
 $(EXAMPLES).o: $(EXAMPLES).cpp $(TARGET_DEPS) $(QLAYER_DEPS) $(EXAMPLES_DEPS)

examples/qAlgorithms.cpp

@@ -5,7 +5,7 @@ QubitLayer grover(unsigned int numQubits, int unsigned dSolution, int unsigned n
 {
     QubitLayer q(numQubits);
     unsigned long long int numStates = q.getNumStates();
-    if (dSolution > numQubits)
+    if (dSolution > numStates)
     {
         std::cout << "\033[31;31m[Error]\033[m" << std::endl;
         std::cout << "Number of states: " << numStates << std::endl;

src/QubitLayer.cpp

@@ -2,6 +2,7 @@
 #include <iostream>
 #include <cmath>
 #include <algorithm>
+#include <omp.h>
 #include "QubitLayer.hpp"
 
 QubitLayer::QubitLayer(unsigned int numQubits, qubitLayer *qL)
@@ -88,8 +89,7 @@ void QubitLayer::applyPauliZ(int target)
 
 void QubitLayer::applyHadamard(int target)
 {
-// map |1> to -hadamardCoef*|1> and |0> to hadamardCoef*|0>
-#pragma omp parallel for shared(qOdd_, qEven_)
+    // map |1> to -hadamardCoef*|1> and |0> to hadamardCoef*|0>
     for (unsigned long long int i = 0; i < numStates; i++)
         if (checkZeroState(i))
         {
@@ -98,25 +98,9 @@ void QubitLayer::applyHadamard(int target)
                 parity ? qOdd_[i] -= hadamardCoef * qEven_[i] : qEven_[i] -= hadamardCoef * qOdd_[i];
             else
                 parity ? qOdd_[i] += hadamardCoef * qEven_[i] : qEven_[i] += hadamardCoef * qOdd_[i];
-            if (!isParallel)
-            {
-                state.flip(target);
-                parity ? qOdd_[state.to_ulong()] += hadamardCoef * qEven_[i] : qEven_[state.to_ulong()] += hadamardCoef * qOdd_[i];
-            }
+            state.flip(target);
+            parity ? qOdd_[state.to_ulong()] += hadamardCoef * qEven_[i] : qEven_[state.to_ulong()] += hadamardCoef * qOdd_[i];
         }
-    if (isParallel)
-    {
-#pragma omp barrier
-// map |0> to hadamardCoef*|1> and |1> to hadamardCoef*|0>
-#pragma omp parallel for shared(qOdd_, qEven_)
-        for (unsigned long long int i = 0; i < numStates; i++)
-            if (checkZeroState(i))
-            {
-                std::bitset<maxQubits> state = i;
-                state.flip(target);
-                parity ? qOdd_[state.to_ulong()] += hadamardCoef * qEven_[i] : qEven_[state.to_ulong()] += hadamardCoef * qOdd_[i];
-            }
-    }
     updateLayer();
 }
 
@@ -125,32 +109,15 @@ void QubitLayer::applyRx(int target, precision theta)
     // compute the sine and cosine of the rotation angle
     precision cosTheta = cos(theta / 2);
     precision sinTheta = sin(theta / 2);
-// map |0> to cosTheta*|0> and |1> to cosTheta*|1>
-#pragma omp parallel for shared(qOdd_, qEven_)
+    // map |0> to cosTheta*|0> and |1> to cosTheta*|1>
     for (unsigned long long int i = 0; i < numStates; i++)
         if (checkZeroState(i))
         {
             parity ? qOdd_[i] += cosTheta * qEven_[i] : qEven_[i] += cosTheta * qOdd_[i];
-            if (!isParallel)
-            {
-                std::bitset<maxQubits> state = i;
-                state.flip(target);
-                parity ? qOdd_[state.to_ulong()] += -complexImg * sinTheta * qEven_[i] : qEven_[state.to_ulong()] += -complexImg * sinTheta * qOdd_[i];
-            }
+            std::bitset<maxQubits> state = i;
+            state.flip(target);
+            parity ? qOdd_[state.to_ulong()] += -complexImg * sinTheta * qEven_[i] : qEven_[state.to_ulong()] += -complexImg * sinTheta * qOdd_[i];
         }
-    if (isParallel)
-    {
-#pragma omp barrier
-// map |1> to -isinTheta*|0> and |0> to -isineTheta*|1>
-#pragma omp parallel for shared(qOdd_, qEven_)
-        for (unsigned long long int i = 0; i < numStates; i++)
-            if (checkZeroState(i))
-            {
-                std::bitset<maxQubits> state = i;
-                state.flip(target);
-                parity ? qOdd_[state.to_ulong()] += -complexImg * sinTheta * qEven_[i] : qEven_[state.to_ulong()] += -complexImg * sinTheta * qOdd_[i];
-            }
-    }
     updateLayer();
 }
 
@@ -159,42 +126,20 @@ void QubitLayer::applyRy(int target, precision theta)
     // compute the sine and cosine of the rotation angle
     precision cosTheta = cos(theta / 2);
     precision sinTheta = sin(theta / 2);
-// map |1> to cosTheta*|1> and |0> to cosTheta*|0>
-#pragma omp parallel for shared(qOdd_, qEven_)
+    // map |1> to cosTheta*|1> and |0> to cosTheta*|0>
     for (unsigned long long int i = 0; i < numStates; i++)
         if (checkZeroState(i))
         {
             parity ? qOdd_[i] += cosTheta * qEven_[i] : qEven_[i] += cosTheta * qOdd_[i];
-            if (!isParallel)
-            {
-                std::bitset<maxQubits> state = i;
-                state.flip(target);
-                // action if bit is 1 (i.e. set)
-                if (state.test(target))
-                    parity ? qOdd_[state.to_ulong()] += sinTheta * qEven_[i] : qEven_[state.to_ulong()] += sinTheta * qOdd_[i];
-                // action if bit is 0 (i.e. not set)
-                else
-                    parity ? qOdd_[state.to_ulong()] -= sinTheta * qEven_[i] : qEven_[state.to_ulong()] -= sinTheta * qOdd_[i];
-            }
+            std::bitset<maxQubits> state = i;
+            state.flip(target);
+            // action if bit is 1 (i.e. set)
+            if (state.test(target))
+                parity ? qOdd_[state.to_ulong()] += sinTheta * qEven_[i] : qEven_[state.to_ulong()] += sinTheta * qOdd_[i];
+            // action if bit is 0 (i.e. not set)
+            else
+                parity ? qOdd_[state.to_ulong()] -= sinTheta * qEven_[i] : qEven_[state.to_ulong()] -= sinTheta * qOdd_[i];
         }
-    if (isParallel)
-    {
-#pragma omp barrier
-// map |0> to sinTheta*|1> and |1> to -sinTheta*|0>
-#pragma omp parallel for shared(qOdd_, qEven_)
-        for (unsigned long long int i = 0; i < numStates; i++)
-            if (checkZeroState(i))
-            {
-                std::bitset<maxQubits> state = i;
-                state.flip(target);
-                // action if bit is 1 (i.e. set)
-                if (state.test(target))
-                    parity ? qOdd_[state.to_ulong()] += sinTheta * qEven_[i] : qEven_[state.to_ulong()] += sinTheta * qOdd_[i];
-                // action if bit is 0 (i.e. not set)
-                else
-                    parity ? qOdd_[state.to_ulong()] -= sinTheta * qEven_[i] : qEven_[state.to_ulong()] -= sinTheta * qOdd_[i];
-            }
-    }
     updateLayer();
 }
 
@@ -203,7 +148,6 @@ void QubitLayer::applyRz(int target, precision theta)
     // compute the sine and cosine of the rotation angle
     precision cosTheta = cos(theta / 2);
     precision sinTheta = sin(theta / 2);
-#pragma omp parallel for shared(qOdd_, qEven_)
     for (unsigned long long int i = 0; i < numStates; i++)
         if (checkZeroState(i))
         {

src/definitions.hpp

@@ -7,8 +7,7 @@ constexpr precision pi{3.14159265358979323846};
 constexpr std::complex<precision> hadamardCoef{0.707106781186548, 0};
 constexpr std::complex<precision> complexImg{0, 1};
 constexpr std::complex<precision> zeroComplex{0, 0};
-constexpr unsigned int maxQubits{32}; //max states allowed by QuantumSim is 2^32
+constexpr unsigned int maxQubits{32}; // max states allowed by QuantumSim is 2^32
 typedef std::complex<precision> qubitLayer;
-static bool isParallel{false};
 
 #endif

src/main.cpp

@@ -6,16 +6,10 @@
 
 int main(int argc, char *argv[])
 {
-    for (int i = 0; i < argc; i++)
-        if (std::string(argv[i]) == "-p")
-        {
-            isParallel = true;
-            std::cout << "Running parallel version of QuantumSim" << std::endl;
-        }
     auto start = std::chrono::steady_clock::now();
-    QubitLayer q = grover(2, 0);
+    QubitLayer q = grover(4, 0);
     auto stop = std::chrono::steady_clock::now();
     auto duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start).count();
     std::cout << "Execution time: " << duration << " µs" << std::endl;
-    q.printMeasurement();
+    q.printQubits();
 }

tests/tests.cpp

@@ -93,14 +93,6 @@ bool testGate(Gates gate)
 
 int main(int argc, char *argv[])
 {
-    for (int i = 0; i < argc; i++)
-        if (std::string(argv[i]) == "-p")
-        {
-            isParallel = true;
-            std::cout << "Testing parallel version of QuantumSim" << std::endl;
-        }
-    if (!isParallel)
-        std::cout << "Testing sequential version of QuantumSim" << std::endl;
     // define variable to store result of the tests
     bool testResult = true;
     std::cout << "\033[34;34m===========Test Results===========\033[m" << std::endl;