added setQuregToWeightedSum (and ToMixture)

Author: TysonRayJones

quest/include/initialisations.h

@@ -116,6 +116,16 @@ void setDensityQuregFlatAmps(Qureg qureg, qindex startInd, qcomp* amps, qindex n
 void setQuregToClone(Qureg targetQureg, Qureg copyQureg);
 
 
+/// @notyetdoced
+/// @notyettested
+void setQuregToWeightedSum(Qureg out, qcomp* coeffs, Qureg* in, int numIn);
+
+
+/// @notyetdoced
+/// @notyettested
+void setQuregToMixture(Qureg out, qreal* probs, Qureg* in, int numIn);
+
+
 /** @notyetdoced
  * @notyettested
  * 
@@ -227,6 +237,20 @@ void setQuregToPartialTrace(Qureg out, Qureg in, std::vector<int> traceOutQubits
 void setQuregToReducedDensityMatrix(Qureg out, Qureg in, std::vector<int> retainQubits);
 
 
+/// @ingroup init_amps
+/// @notyetdoced
+/// @cpponly
+/// @see setQuregToWeightedSum()
+void setQuregToWeightedSum(Qureg out, std::vector<qcomp> coeffs, std::vector<Qureg> in);
+
+
+/// @ingroup init_amps
+/// @notyetdoced
+/// @cpponly
+/// @see setQuregToMixture()
+void setQuregToMixture(Qureg out, std::vector<qreal> probs, std::vector<Qureg> in);
+
+
 #endif // __cplusplus
 
 

quest/src/api/decoherence.cpp

@@ -134,7 +134,7 @@ void mixQureg(Qureg outQureg, Qureg inQureg, qreal inProb) {
     validate_quregFields(outQureg, __func__);
     validate_quregFields(inQureg, __func__);
     validate_probability(inProb, __func__);
-    validate_quregsCanBeMixed(outQureg, inQureg, __func__); // checks outQureg is densmatr
+    validate_quregPairCanBeMixed(outQureg, inQureg, __func__); // checks outQureg is densmatr
 
     qreal outProb = 1 - inProb;
     localiser_densmatr_mixQureg(outProb, outQureg, inProb, inQureg);

quest/src/api/initialisations.cpp

@@ -252,6 +252,34 @@ void setQuregToReducedDensityMatrix(Qureg out, Qureg in, int* retainQubits, int
 }
 
 
+void setQuregToWeightedSum(Qureg out, qcomp* coeffs, Qureg* in, int numIn) {
+    validate_quregFields(out, __func__);
+    validate_numQuregsInSum(numIn, __func__);
+    validate_quregsCanBeSummed(out, in, numIn, __func__); // also validates all init
+
+    auto coeffVec = util_getVector(coeffs, numIn);
+    auto inVec = util_getVector(in, numIn);
+    localiser_statevec_setQuregToWeightedSum(out, coeffVec, inVec);
+}
+
+
+void setQuregToMixture(Qureg out, qreal* probs, Qureg* in, int numIn) {
+    validate_quregFields(out, __func__);
+    validate_quregIsDensityMatrix(out, __func__);
+    validate_numQuregsInSum(numIn, __func__);
+    validate_quregsCanBeMixed(out, in, numIn, __func__); // also validates all init & densmatr
+    validate_probabilities(probs, numIn, __func__);
+
+    // convert probs to complex (assume this alloc never fails)
+    vector<qcomp> coeffVec(numIn);
+    for (int i=0; i<numIn; i++)
+        coeffVec[i] = getQcomp(probs[i], 0);
+
+    auto inVec = util_getVector(in, numIn);
+    localiser_statevec_setQuregToWeightedSum(out, coeffVec, inVec);
+}
+
+
 } // end de-mangler
 
 
@@ -296,3 +324,15 @@ void setQuregToPartialTrace(Qureg out, Qureg in, vector<int> traceOutQubits) {
 void setQuregToReducedDensityMatrix(Qureg out, Qureg in, vector<int> retainQubits) {
     setQuregToReducedDensityMatrix(out, in, retainQubits.data(), retainQubits.size());
 }
+
+void setQuregToWeightedSum(Qureg out, vector<qcomp> coeffs, vector<Qureg> in) {
+    validate_numQuregsMatchesCoeffs(in.size(), coeffs.size(), __func__);
+
+    setQuregToWeightedSum(out, coeffs.data(), in.data(), in.size());
+}
+
+void setQuregToMixture(Qureg out, vector<qreal> probs, vector<Qureg> in) {
+    validate_numQuregsMatchesProbs(in.size(), probs.size(), __func__);
+
+    setQuregToMixture(out, probs.data(), in.data(), in.size());
+}

quest/src/core/accelerator.cpp

@@ -74,6 +74,10 @@ using std::min;
 #endif
 
 
+#define GET_FUNC_OPTIMISED_FOR_NUM_QUREGS(f, numquregs) \
+    (vector <decltype(&f<0>)> {&f<0>, &f<1>, &f<2>, &f<3>, &f<4>, &f<5>, &f<-1>}) \
+    [std::min((int) numquregs, MAX_OPTIMISED_NUM_QUREGS + 1)]
+
 #define GET_FUNC_OPTIMISED_FOR_NUM_CTRLS(f, numctrls) \
     (vector <decltype(&f<0>)> {&f<0>, &f<1>, &f<2>, &f<3>, &f<4>, &f<5>, &f<-1>}) \
     [std::min((int) numctrls, MAX_OPTIMISED_NUM_CTRLS + 1)]
@@ -97,6 +101,11 @@ using std::min;
 #define ARR(f) vector<decltype(&f<0,0>)>
 
 
+#define GET_CPU_OR_GPU_FUNC_OPTIMISED_FOR_NUM_QUREGS(funcsuffix, qureg, numquregs) \
+    ((qureg.isGpuAccelerated)? \
+        GET_FUNC_OPTIMISED_FOR_NUM_QUREGS( gpu_##funcsuffix, numquregs ) : \
+        GET_FUNC_OPTIMISED_FOR_NUM_QUREGS( cpu_##funcsuffix, numquregs ))
+
 #define GET_CPU_OR_GPU_FUNC_OPTIMISED_FOR_NUM_CTRLS(funcsuffix, qureg, numctrls) \
     ((qureg.isGpuAccelerated)? \
         GET_FUNC_OPTIMISED_FOR_NUM_CTRLS( gpu_##funcsuffix, numctrls ) : \
@@ -539,6 +548,14 @@ void accel_statevector_anyCtrlAnyTargZOrPhaseGadget_sub(Qureg qureg, vector<int>
  */
 
 
+void accel_statevec_setQuregToWeightedSum_sub(Qureg outQureg, vector<qcomp> coeffs, vector<Qureg> inQuregs) {
+
+    // consult outQureg's deployment since others are prior validated to match
+    auto func = GET_CPU_OR_GPU_FUNC_OPTIMISED_FOR_NUM_QUREGS( statevec_setQuregToWeightedSum_sub, outQureg, inQuregs.size() );
+    func(outQureg, coeffs, inQuregs);
+}
+
+
 void accel_statevec_setQuregToSuperposition_sub(qcomp facOut, Qureg outQureg, qcomp fac1, Qureg inQureg1, qcomp fac2, Qureg inQureg2) {
 
     // consult outQureg's deployment (other quregs should match, though we dangerously do not assert this post-validation)

quest/src/core/accelerator.hpp

@@ -35,7 +35,7 @@ using std::vector;
  * used by cpu_subroutines.cpp and gpu_subroutines to force the compiler
  * to instantiate and compile their template definitions with the given
  * explicit parameters below. Notice the final parameter is always -1, 
- * to handle when the number of controls or targets is not known at 
+ * to handle when the number of parameters (e.g. controls) is not known at 
  * compile-time (it is larger than a bespoke, optimised instantiations), 
  * causing the optimised function to fallback to a suboptimal but general 
  * implementation.
@@ -44,6 +44,7 @@ using std::vector;
 // must match the macros below, and those in accelerator.cpp
 #define MAX_OPTIMISED_NUM_CTRLS 5
 #define MAX_OPTIMISED_NUM_TARGS 5
+#define MAX_OPTIMISED_NUM_QUREGS 5
 
 
 #define INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_TARGS(returntype, funcname, args) \
@@ -58,6 +59,9 @@ using std::vector;
 #define INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_CTRLS(returntype, funcname, args) \
     INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_TARGS(returntype, funcname, args)
 
+#define INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_QUREGS(returntype, funcname, args) \
+    INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_TARGS(returntype, funcname, args)
+
 
 #define INSTANTIATE_FUNC_OPTIMISED_FOR_NUM_CTRLS_AND_TARGS(returntype, funcname, args) \
     private_INSTANTIATE(returntype, funcname, 0, args) \
@@ -228,6 +232,8 @@ void accel_statevector_anyCtrlPauliTensorOrGadget_subB(Qureg qureg, vector<int>
  * QUREG COMBINATION
  */
 
+void accel_statevec_setQuregToWeightedSum_sub(Qureg outQureg, vector<qcomp> coeffs, vector<Qureg> inQuregs);
+
 void accel_statevec_setQuregToSuperposition_sub(qcomp facOut, Qureg outQureg, qcomp fac1, Qureg inQureg1, qcomp fac2, Qureg inQureg2);
 
 void accel_densmatr_mixQureg_subA(qreal outProb, Qureg out, qreal inProb, Qureg in);

quest/src/core/localiser.cpp

@@ -1371,6 +1371,19 @@ void localiser_statevec_anyCtrlPauliGadget(Qureg qureg, vector<int> ctrls, vecto
  */
 
 
+void localiser_statevec_setQuregToWeightedSum(Qureg outQureg, vector<qcomp> coeffs, vector<Qureg> inQuregs) {
+
+    /// @todo
+    /// this function requires (as validated) distributions are identical.
+    /// It would be trivial to generalise this so that Qureg distributions
+    /// can differ (we merely spoof local Quregs, offsetting their memory).
+    /// They must still however be identically GPU-accelerated; this is a
+    /// low priority because this situation is non-sensical
+
+    accel_statevec_setQuregToWeightedSum_sub(outQureg, coeffs, inQuregs);
+}
+
+
 void localiser_statevec_setQuregToSuperposition(qcomp facOut, Qureg outQureg, qcomp fac1, Qureg inQureg1, qcomp fac2, Qureg inQureg2) {
 
     /// @todo

quest/src/core/localiser.hpp

@@ -127,6 +127,8 @@ void localiser_statevec_anyCtrlPhaseGadget(Qureg qureg, vector<int> ctrls, vecto
  * QUREG COMBINATION
  */
 
+void localiser_statevec_setQuregToWeightedSum(Qureg outQureg, vector<qcomp> coeffs, vector<Qureg> inQuregs);
+
 void localiser_statevec_setQuregToSuperposition(qcomp facOut, Qureg outQureg, qcomp fac1, Qureg inQureg1, qcomp fac2, Qureg inQureg2);
 
 void localiser_densmatr_mixQureg(qreal outProb, Qureg out, qreal inProb, Qureg in);

quest/src/core/utilities.cpp

@@ -218,15 +218,6 @@ qindex util_getBitMask(vector<int> ctrls, vector<int> ctrlStates, vector<int> ta
     return util_getBitMask(qubits, states);
 }
 
-vector<int> util_getVector(int* qubits, int numQubits) {
-
-    // permit qubits=nullptr, overriding numQubits (might be non-zero)
-    if (qubits == nullptr)
-        return {};
-
-    return vector<int> (qubits, qubits + numQubits);
-}
-
 
 
 /*
@@ -384,7 +375,7 @@ qreal util_getSum(vector<qreal> list) {
     qreal sum = 0;
     qreal y, t, c=0;
 
-    // complex Kahan summation
+    // Kahan summation
     for (auto& x : list) {
         y = x - c;
         t = sum + y;
@@ -1118,6 +1109,23 @@ qreal util_getMaxProbOfTwoQubitDepolarising() {
  * TEMPORARY MEMORY ALLOCATION
  */
 
+template <typename T>
+vector<T> getVector(T* ptr, int length) {
+
+    // permit nullptr to indicate empty list, regardless of length
+    if (ptr == nullptr)
+        return {};
+
+    // assumes memory alloc failure is impossible
+    return vector<T> (ptr, ptr + length);
+}
+
+vector<int>   util_getVector(int*   ptr, int length) { return getVector(ptr, length); }
+vector<qreal> util_getVector(qreal* ptr, int length) { return getVector(ptr, length); }
+vector<qcomp> util_getVector(qcomp* ptr, int length) { return getVector(ptr, length); }
+vector<Qureg> util_getVector(Qureg* ptr, int length) { return getVector(ptr, length); }
+
+
 template <typename T>
 void tryAllocVector(vector<T> &vec, qindex size, std::function<void()> errFunc) {
 

quest/src/core/utilities.hpp

@@ -59,8 +59,6 @@ vector<int> util_getBraQubits(vector<int> ketQubits, Qureg qureg);
 
 vector<int> util_getNonTargetedQubits(int* targets, int numTargets, int numQubits);
 
-vector<int> util_getVector(int* qubits, int numQubits);
-
 vector<int> util_getConcatenated(vector<int> list1, vector<int> list2);
 
 vector<int> util_getSorted(vector<int> list);
@@ -401,6 +399,13 @@ qreal util_getMaxProbOfTwoQubitDepolarising();
  * TEMPORARY MEMORY ALLOCATION
  */
 
+// alloc assumed to never fail
+vector<int>   util_getVector(int*   ptr, int length);
+vector<qreal> util_getVector(qreal* ptr, int length);
+vector<qcomp> util_getVector(qcomp* ptr, int length);
+vector<Qureg> util_getVector(Qureg* ptr, int length);
+
+// calls errFunc when alloc fails
 void util_tryAllocVector(vector<qreal>    &vec, qindex size, std::function<void()> errFunc);
 void util_tryAllocVector(vector<qcomp>    &vec, qindex size, std::function<void()> errFunc);
 void util_tryAllocVector(vector<qcomp*>   &vec, qindex size, std::function<void()> errFunc);