scaffolded HPC deployment

Author: TysonRayJones

quest/include/operations.h

@@ -7,6 +7,27 @@
 #ifndef OPERATIONS_H
 #define OPERATIONS_H
 
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
 
+// enable invocation by both C and C++ binaries
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+
+
+/*
+ * GATES
+ */
+
+void unitary(Qureg qureg, int target, CompMatr1 matrix);
+
+
+
+// end de-mangler
+#ifdef __cplusplus
+}
+#endif
 
 #endif // OPERATIONS_H

quest/src/api/operations.cpp

@@ -2,4 +2,28 @@
  * API definitions for effecting operators (such as gates and unitaries) 
  * upon Quregs which are instantiated as either statevectors or 
  * density matrices.
- */
+ */
+
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+#include "quest/src/core/validation.hpp"
+#include "quest/src/core/utilities.hpp"
+#include "quest/src/core/localiser.hpp"
+
+
+// enable invocation by both C and C++ binaries
+extern "C" {
+
+
+void unitary(Qureg qureg, int target, CompMatr1 matrix) {
+    validate_target(qureg, target, __func__);
+    validate_matrixIsUnitary(matrix, __func__);
+
+    statevec_oneTargetGate(qureg, target, matrix);
+    if (qureg.isDensityMatrix)
+        statevec_oneTargetGate(qureg, util_getShifted(target, qureg), util_getConj(matrix));
+}
+
+
+} // end de-mangler

quest/src/core/accelerator.cpp

@@ -6,4 +6,28 @@
  * and location (qureg's amplitudes or buffer space), as is
  * performed by localiser.cpp. These subroutines are ergo
  * embarrassingly parallel.
- */
+ */
+
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+#include "quest/src/cpu/cpu_subroutines.hpp"
+#include "quest/src/gpu/gpu_subroutines.hpp"
+
+
+
+void statevec_oneTargetGate_subA(Qureg qureg, int target, CompMatr1 matrix) {
+
+    if (qureg.isGpuAccelerated)
+        gpu_statevec_oneTargetGate_subA(qureg, target, matrix);
+    else
+        cpu_statevec_oneTargetGate_subA(qureg, target, matrix);
+
+}
+void statevec_oneTargetGate_subB(Qureg qureg, qcomp fac0, qcomp fac1) {
+
+    if (qureg.isGpuAccelerated)
+        gpu_statevec_oneTargetGate_subB(qureg, fac0, fac1);
+    else
+        cpu_statevec_oneTargetGate_subB(qureg, fac0, fac1);
+}

quest/src/core/accelerator.hpp

@@ -11,6 +11,15 @@
 #ifndef ACCELERATOR_HPP
 #define ACCELERATOR_HPP
 
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+
+
+void statevec_oneTargetGate_subA(Qureg qureg, int target, CompMatr1 matrix);
+
+void statevec_oneTargetGate_subB(Qureg qureg, qcomp fac0, qcomp fac1);
+
 
 
 #endif // ACCELERATOR_HPP

quest/src/core/localiser.cpp

@@ -6,4 +6,53 @@
  * embarrassingly parallel subroutines in accelerator.cpp. This is
  * done agnostically of whether amplitudes of the Qureg are being
  * stored in RAM (CPU) or VRAM (GPU).
- */
+ */
+
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+#include "quest/src/comm/comm_routines.hpp"
+#include "quest/src/core/bitwise.hpp"
+#include "quest/src/core/accelerator.hpp"
+
+
+
+/*
+ * PRIVATE
+ */
+
+
+bool oneTargetGateRequiresComm(Qureg qureg, int target) {
+    if (!qureg.isDistributed)
+        return false;
+    
+    return target >= qureg.logNumAmpsPerNode;
+}
+
+
+
+/*
+ * OPERATORS
+ */
+
+
+void statevec_oneTargetGate(Qureg qureg, int target, CompMatr1 matrix) {
+
+    if (!oneTargetGateRequiresComm(qureg, target))
+        statevec_oneTargetGate_subA(qureg, target, matrix);
+
+    else {
+        // exchange all amps (receive to buffer)
+        int rankTarget = target - qureg.logNumAmpsPerNode;
+        int pairRank = flipBit(qureg.rank, rankTarget);
+
+        comm_exchangeAmpsToBuffers(qureg, pairRank);
+
+        // extract relevant gate elements
+        int bit = getBit(qureg.rank, rankTarget);
+        qcomp fac0 = matrix.elems[bit][bit];
+        qcomp fac1 = matrix.elems[bit][!bit];
+
+        statevec_oneTargetGate_subB(qureg, fac0, fac1);
+    }
+}

quest/src/core/localiser.hpp

@@ -11,6 +11,17 @@
 #ifndef LOCALISER_HPP
 #define LOCALISER_HPP
 
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+
+
+/*
+ * OPERATORS
+ */
+
+void statevec_oneTargetGate(Qureg qureg, int target, CompMatr1 matrix);
+
 
 
 #endif // LOCALISER_HPP

quest/src/core/validation.cpp

@@ -11,6 +11,7 @@
 #include "quest/src/core/errors.hpp"
 #include "quest/src/core/bitwise.hpp"
 #include "quest/src/core/memory.hpp"
+#include "quest/src/core/utilities.hpp"
 #include "quest/src/comm/comm_config.hpp"
 #include "quest/src/cpu/cpu_config.hpp"
 #include "quest/src/gpu/gpu_config.hpp"
@@ -180,8 +181,17 @@ namespace report {
     std::string INVALID_EXISTING_ALLOC_OF_COMPLEX_MATRIX =
         "Invalid CompMatrN. One or more rows of the 2D elements array was seemingly unallocated.";
 
-    std::string INVALID_COMP_MATR_PROPERTIES = 
-        "Invalid CompMatrN. Targeted ${NUM_QUBITS} qubits and had a dimension of ${NUM_ROWS} x ${NUM_ROWS}.";
+    std::string INVALID_COMP_MATR_1_FIELDS =
+        "Invalid CompMatr1. Targeted ${NUM_QUBITS} qubits (instead of 1) and had a dimension of ${NUM_ROWS}x${NUM_ROWS} (instead of 2x2). It is likely this matrix was not initialised with getCompMatr1().";
+
+    std::string INVALID_COMP_MATR_2_FIELDS =
+        "Invalid CompMatr2. Targeted ${NUM_QUBITS} qubits (instead of 2) and had a dimension of ${NUM_ROWS}x${NUM_ROWS} (instead of 4x4). It is likely this matrix was not initialised with getCompMatr2().";
+
+    std::string INVALID_COMP_MATR_N_FIELDS =
+        "Invalid CompMatrN. Targeted ${NUM_QUBITS} qubits and had a dimension of ${NUM_ROWS}x${NUM_ROWS}. It is likely this matrix was not initialised with createCompMatrN().";
+
+    std::string MATRIX_NOT_UNITARY = 
+        "The given complex matrix was not (approximately) unitary.";
 
 
     /*
@@ -190,6 +200,14 @@ namespace report {
 
     std::string INVALID_STATE_INDEX = 
         "Classical state index ${STATE_IND} is invalid for the given ${NUM_QUBITS} qubit Qureg. Index must be greater than or equal to zero, and cannot equal nor exceed the number of unique classical states (2^${NUM_QUBITS} = ${NUM_STATES}).";
+
+    
+    /*
+     * OPERATOR PARAMETERS
+     */
+    
+    std::string INVALID_TARGET_QUBIT = 
+        "Invalid target qubit (${TARGET}). Must be greater than or equal to zero, and less than the number of qubits in the Qureg (${NUM_QUBITS}).";
 }
 
 
@@ -204,7 +222,6 @@ void default_invalidQuESTInputError(const char* msg, const char* func) {
     // safe to call even before MPI has been setup
     if (comm_getRank() == 0)
         std::cout 
-            << std::endl
             << "QuEST encountered a validation error during function '" << func << "':\n"
             << msg << "\nExiting..." 
             << std::endl;
@@ -672,20 +689,52 @@ void validate_newOrExistingMatrixAllocs(CompMatrN matr, bool isNewMatr, const ch
         assertThat(matr.elems[r] != NULL, errMsg, caller);
 }
 
-void validate_matrixInit(CompMatrN matr, const char* caller) {
+
+/*
+ * EXISTING MATRICES
+ */
+
+template <class T>
+void assertMatrixFieldsAreValid(T matr, std::string errMsg, const char* caller) {
 
     tokenSubs vars = {
         {"${NUM_QUBITS}", matr.numQubits},
         {"${NUM_ROWS}",   matr.numRows}};
 
     bool validDims = (matr.numQubits >= 1) && (matr.numRows == powerOf2(matr.numQubits));
-    assertThat(validDims, report::INVALID_COMP_MATR_PROPERTIES, vars, caller);
+    assertThat(validDims, errMsg, vars, caller);
+}
+
+void validate_matrixInit(CompMatr1 matr, const char* caller) {
+    assertMatrixFieldsAreValid(matr, report::INVALID_COMP_MATR_1_FIELDS, caller);
+}
+void validate_matrixInit(CompMatr2 matr, const char* caller) {
+    assertMatrixFieldsAreValid(matr, report::INVALID_COMP_MATR_2_FIELDS, caller);
+}
+void validate_matrixInit(CompMatrN matr, const char* caller) {
+    assertMatrixFieldsAreValid(matr, report::INVALID_COMP_MATR_N_FIELDS, caller);
+    validate_newOrExistingMatrixAllocs(matr, false, caller);
+}
+
+template <class T> 
+void assertMatrixIsUnitary(T matr, const char* caller) {
 
-    bool isNewMatr = false;
-    validate_newOrExistingMatrixAllocs(matr, isNewMatr, caller);
+    validate_matrixInit(matr, caller);
+    assertThat(util_isUnitary(matr), report::MATRIX_NOT_UNITARY, caller);
+}
+
+void validate_matrixIsUnitary(CompMatr1 matr, const char* caller) {
+    assertMatrixIsUnitary(matr, caller);
+}
+void validate_matrixIsUnitary(CompMatr2 matr, const char* caller) {
+    assertMatrixIsUnitary(matr, caller);
+}
+void validate_matrixIsUnitary(CompMatrN matr, const char* caller) {
+    assertMatrixIsUnitary(matr, caller);
 }
 
 
+
 /*
  * QUREG INITIALISATIONS
  */
@@ -701,3 +750,18 @@ void validate_initClassicalStateIndex(Qureg qureg, qindex ind, const char* calle
 
     assertThat(ind >= 0 && ind < maxIndExcl, report::INVALID_STATE_INDEX, vars, caller);
 }
+
+
+
+/*
+ * OPERATOR PARAMETERS
+ */
+
+void validate_target(Qureg qureg, int target, const char* caller) {
+
+    tokenSubs vars = {
+        {"${TARGET}",  target},
+        {"${NUM_QUBITS}", qureg.numQubits}};
+
+    assertThat(target >= 0 && target < qureg.numQubits, report::INVALID_TARGET_QUBIT, vars, caller);
+}

quest/src/core/validation.hpp

@@ -66,8 +66,14 @@ void validate_newOrExistingMatrixAllocs(CompMatrN matr, bool isNewMatr, const ch
  * EXISTING MATRIX
  */
 
+void validate_matrixInit(CompMatr1 matr, const char* caller);
+void validate_matrixInit(CompMatr2 matr, const char* caller);
 void validate_matrixInit(CompMatrN matr, const char* caller);
 
+void validate_matrixIsUnitary(CompMatr1 matr, const char* caller);
+void validate_matrixIsUnitary(CompMatr2 matr, const char* caller);
+void validate_matrixIsUnitary(CompMatrN matr, const char* caller);
+
 
 
 /*
@@ -78,4 +84,12 @@ void validate_initClassicalStateIndex(Qureg qureg, qindex ind, const char* calle
 
 
 
+/*
+ * OPERATOR PARAMETERS
+ */
+
+void validate_target(Qureg qureg, int target, const char* caller);
+
+
+
 #endif // VALIDATION_HPP

quest/src/cpu/cpu_subroutines.cpp

@@ -3,4 +3,40 @@
  * accelerator.cpp. 
  */
 
-#include "quest/src/cpu/cpu_config.hpp"
+#include "quest/include/types.h"
+#include "quest/include/qureg.h"
+#include "quest/include/structures.h"
+
+#include "quest/src/core/bitwise.hpp"
+
+
+
+void cpu_statevec_oneTargetGate_subA(Qureg qureg, int target, CompMatr1 matrix) {
+
+    qindex numIts = qureg.numAmpsPerNode / 2;
+
+    qcomp m00 = matrix.elems[0][0];
+    qcomp m01 = matrix.elems[0][1];
+    qcomp m10 = matrix.elems[1][0];
+    qcomp m11 = matrix.elems[1][1];
+    
+    for (qindex j=0; j<numIts; j++) {
+
+        qindex i0 = insertBit(j, target, 0);
+        qindex i1 = flipBit(i0, target);
+        
+        qcomp amp0 = qureg.cpuAmps[i0];
+        qcomp amp1 = qureg.cpuAmps[i1];
+        
+        qureg.cpuAmps[i0] = m00*amp0 + m01*amp1;
+        qureg.cpuAmps[i1] = m10*amp0 + m11*amp1;
+    }
+}
+
+void cpu_statevec_oneTargetGate_subB(Qureg qureg, qcomp fac0, qcomp fac1) {
+
+    qindex numIts = qureg.numAmpsPerNode;
+
+    for (qindex i=0; i<numIts; i++)
+        qureg.cpuAmps[i] = fac0*qureg.cpuAmps[i] + fac1*qureg.cpuCommBuffer[i];
+}