operations.cpp

/** @file
 * API definitions for effecting operators (including
 * unitaries, projectors, channels, Hermitians, and
 * arbitrary matrices) upon Quregs, which can be
 * statevectors or density matrices.
 * 
 * @author Tyson Jones
 */
 
#include "quest/include/qureg.h"
#include "quest/include/paulis.h"
#include "quest/include/matrices.h"
#include "quest/include/operations.h"
#include "quest/include/calculations.h"
 
#include "quest/src/core/validation.hpp"
#include "quest/src/core/utilities.hpp"
#include "quest/src/core/randomiser.hpp"
#include "quest/src/core/localiser.hpp"
#include "quest/src/core/bitwise.hpp"
#include "quest/src/core/constants.hpp"
#include "quest/src/core/paulilogic.hpp"
 
#include <vector>
 
using std::vector;
 
 
 
/*
 * PRVIATE UTILITIES
 */
 
// T can be CompMatr, CompMatr1, CompMatr2, DiagMatr, DiagMatr1, DiagMatr2
template <class T>
void validateAndApplyAnyCtrlAnyTargUnitaryMatrix(Qureg qureg, int* ctrls, int* states, int numCtrls, int* targs, int numTargs, T matr, const char* caller) {
    validate_quregFields(qureg, caller);
    validate_controlsAndTargets(qureg, ctrls, numCtrls, targs, numTargs, caller);
    validate_controlStates(states, numCtrls, caller);
    validate_matrixDimMatchesTargets(matr, numTargs, caller); // also checks fields and is-synced
    validate_matrixIsUnitary(matr, caller); // harmlessly rechecks fields and is-synced
    if (util_isDenseMatrixType<T>())
        validate_mixedAmpsFitInNode(qureg, numTargs, caller);
 
    List64 ctrlList  = lists_getList64(ctrls,  numCtrls);
    List64 stateList = util_getList64OrAllOnes(states, numCtrls);
    List64 targList  = lists_getList64(targs,  numTargs);
 
    bool conj = false;
    localiser_statevec_anyCtrlAnyTargAnyMatr(qureg, ctrlList, stateList, targList, matr, conj);
 
    if (!qureg.isDensityMatrix)
        return;
 
    conj = true;
    ctrlList = util_getBraQubits(ctrlList, qureg);
    targList = util_getBraQubits(targList, qureg);
    localiser_statevec_anyCtrlAnyTargAnyMatr(qureg, ctrlList, stateList, targList, matr, conj);
 
    /// @todo
    /// the above logic always performs two in-turn operations upon density matrices, 
    /// though when matr is diagonal (DiagMatr*), they can be trivially combined into 
    /// a single operation which enumerates the state only once. We perform this
    /// optimisation for FullStateDiagMatr elsewhere. Consider optimising here too!
}
 
 
 
/*
 * CompMatr1
 */
 
extern "C" {
 
void applyCompMatr1(Qureg qureg, int target, CompMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, &target, 1, matrix, __func__);
}
 
void applyControlledCompMatr1(Qureg qureg, int control, int target, CompMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, &target, 1, matrix, __func__);
}
 
void applyMultiControlledCompMatr1(Qureg qureg, int* controls, int numControls, int target, CompMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledCompMatr1(Qureg qureg, int* controls, int* states, int numControls, int target, CompMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledCompMatr1(Qureg qureg, vector<int> controls, int target, CompMatr1 matrix) {
 
    applyMultiControlledCompMatr1(qureg, controls.data(), controls.size(), target, matrix);
}
 
void applyMultiStateControlledCompMatr1(Qureg qureg, vector<int> controls, vector<int> states, int target, CompMatr1 matrix) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledCompMatr1(qureg, controls.data(), states.data(), controls.size(), target, matrix);
}
 
 
 
/*
 * CompMatr2
 */
 
extern "C" {
 
void applyCompMatr2(Qureg qureg, int target1, int target2, CompMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, targs, 2, matrix, __func__);
}
 
void applyControlledCompMatr2(Qureg qureg, int control, int target1, int target2, CompMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, targs, 2, matrix, __func__);
}
 
void applyMultiControlledCompMatr2(Qureg qureg, int* controls, int numControls, int target1, int target2, CompMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, targs, 2, matrix, __func__);
}
 
void applyMultiStateControlledCompMatr2(Qureg qureg, int* controls, int* states, int numControls, int target1, int target2, CompMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, targs, 2, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledCompMatr2(Qureg qureg, vector<int> controls, int target1, int target2, CompMatr2 matr) {
 
    applyMultiControlledCompMatr2(qureg, controls.data(), controls.size(), target1, target2, matr);
}
 
void applyMultiStateControlledCompMatr2(Qureg qureg, vector<int> controls, vector<int> states, int target1, int target2, CompMatr2 matr) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledCompMatr2(qureg, controls.data(), states.data(), controls.size(), target1, target2, matr);
}
 
 
 
/*
 * CompMatr
 */
 
extern "C" {
 
void applyCompMatr(Qureg qureg, int* targets, int numTargets, CompMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, targets, numTargets, matrix, __func__);
}
 
void applyControlledCompMatr(Qureg qureg, int control, int* targets, int numTargets, CompMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, targets, numTargets, matrix, __func__);
}
 
void applyMultiControlledCompMatr(Qureg qureg, int* controls, int numControls, int* targets, int numTargets, CompMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, targets, numTargets, matrix, __func__);
}
 
void applyMultiStateControlledCompMatr(Qureg qureg, int* controls, int* states, int numControls, int* targets, int numTargets, CompMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, targets, numTargets, matrix, __func__);
}
 
} // end de-mangler
 
void applyCompMatr(Qureg qureg, vector<int> targets, CompMatr matr) {
 
    applyCompMatr(qureg, targets.data(), targets.size(), matr);
}
 
void applyControlledCompMatr(Qureg qureg, int control, vector<int> targets, CompMatr matr) {
 
    applyControlledCompMatr(qureg, control, targets.data(), targets.size(), matr);
}
 
void applyMultiControlledCompMatr(Qureg qureg, vector<int> controls, vector<int> targets, CompMatr matr) {
 
    applyMultiControlledCompMatr(qureg, controls.data(), controls.size(), targets.data(), targets.size(), matr);
}
 
void applyMultiStateControlledCompMatr(Qureg qureg, vector<int> controls, vector<int> states, vector<int> targets, CompMatr matr) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledCompMatr(qureg, controls.data(), states.data(), controls.size(), targets.data(), targets.size(), matr);
}
 
 
 
/*
 * DiagMatr1
 */
 
extern "C" {
 
void applyDiagMatr1(Qureg qureg, int target, DiagMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, &target, 1, matrix, __func__);
}
 
void applyControlledDiagMatr1(Qureg qureg, int control, int target, DiagMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, &target, 1, matrix, __func__);
}
 
void applyMultiControlledDiagMatr1(Qureg qureg, int* controls, int numControls, int target, DiagMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledDiagMatr1(Qureg qureg, int* controls, int* states, int numControls, int target, DiagMatr1 matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledDiagMatr1(Qureg qureg, vector<int> controls, int target, DiagMatr1 matr) {
 
    applyMultiControlledDiagMatr1(qureg, controls.data(), controls.size(), target, matr);
}
 
void applyMultiStateControlledDiagMatr1(Qureg qureg, vector<int> controls, vector<int> states, int target, DiagMatr1 matr) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledDiagMatr1(qureg, controls.data(), states.data(), controls.size(), target, matr);
}
 
 
 
/*
 * DiagMatr2
 */
 
extern "C" {
 
void applyDiagMatr2(Qureg qureg, int target1, int target2, DiagMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, targs, 2, matrix, __func__);
}
 
void applyControlledDiagMatr2(Qureg qureg, int control, int target1, int target2, DiagMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, targs, 2, matrix, __func__);
}
 
void applyMultiControlledDiagMatr2(Qureg qureg, int* controls, int numControls, int target1, int target2, DiagMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, targs, 2, matrix, __func__);
}
 
void applyMultiStateControlledDiagMatr2(Qureg qureg, int* controls, int* states, int numControls, int target1, int target2, DiagMatr2 matrix) {
 
    int targs[] = {target1, target2};
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, targs, 2, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledDiagMatr2(Qureg qureg, vector<int> controls, int target1, int target2, DiagMatr2 matr) {
 
    applyMultiControlledDiagMatr2(qureg, controls.data(), controls.size(), target1, target2, matr);
}
 
void applyMultiStateControlledDiagMatr2(Qureg qureg, vector<int> controls, vector<int> states, int target1, int target2, DiagMatr2 matr) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledDiagMatr2(qureg, controls.data(), states.data(), controls.size(), target1, target2, matr);
}
 
 
 
/*
 * DiagMatr
 */
 
extern "C" {
 
void applyDiagMatr(Qureg qureg, int* targets, int numTargets, DiagMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, nullptr, nullptr, 0, targets, numTargets, matrix, __func__);
}
 
void applyControlledDiagMatr(Qureg qureg, int control, int* targets, int numTargets, DiagMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, &control, nullptr, 1, targets, numTargets, matrix, __func__);
}
 
void applyMultiControlledDiagMatr(Qureg qureg, int* controls, int numControls, int* targets, int numTargets, DiagMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, nullptr, numControls, targets, numTargets, matrix, __func__);
}
 
void applyMultiStateControlledDiagMatr(Qureg qureg, int* controls, int* states, int numControls, int* targets, int numTargets, DiagMatr matrix) {
 
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, targets, numTargets, matrix, __func__);
}
 
} // end de-mangler
 
void applyDiagMatr(Qureg qureg, vector<int> targets, DiagMatr matrix) {
 
    applyDiagMatr(qureg, targets.data(), targets.size(), matrix);
}
 
void applyControlledDiagMatr(Qureg qureg, int control, vector<int> targets, DiagMatr matrix) {
 
    applyControlledDiagMatr(qureg, control, targets.data(), targets.size(), matrix);
}
 
void applyMultiControlledDiagMatr(Qureg qureg, vector<int> controls, vector<int> targets, DiagMatr matrix) {
 
    applyMultiControlledDiagMatr(qureg, controls.data(), controls.size(), targets.data(), targets.size(), matrix);
}
 
void applyMultiStateControlledDiagMatr(Qureg qureg, vector<int> controls, vector<int> states, vector<int> targets, DiagMatr matrix) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
    
    applyMultiStateControlledDiagMatr(qureg, controls.data(), states.data(), controls.size(), targets.data(), targets.size(), matrix);
}
 
 
 
/*
 * DiagMatrPower
 *
 * which still assert unitarity despite that passing
 * a non-real exponent is permitted
 */
 
extern "C" {
 
void applyDiagMatrPower(Qureg qureg, int* targets, int numTargets, DiagMatr matrix, qcomp exponent)  {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
    validate_matrixDimMatchesTargets(matrix, numTargets, __func__); // also checks fields and is-synced
    validate_matrixIsUnitary(matrix, __func__);
    validate_matrixExpIsNonDiverging(matrix, exponent, __func__);
    validate_unitaryExponentIsReal(exponent, __func__);
 
    // notice exponent is a 'qcomp' not a 'qreal' despite validation; see applyMultiStateControlledDiagMatrPower()
 
    // harmlessly re-validates
    applyMultiStateControlledDiagMatrPower(qureg, nullptr, nullptr, 0, targets, numTargets, matrix, exponent);
}
 
void applyControlledDiagMatrPower(Qureg qureg, int control, int* targets, int numTargets, DiagMatr matrix, qcomp exponent) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTargets(qureg, control, targets, numTargets, __func__);
    validate_matrixDimMatchesTargets(matrix, numTargets, __func__); // also checks fields and is-synced
    validate_matrixIsUnitary(matrix, __func__);
    validate_matrixExpIsNonDiverging(matrix, exponent, __func__);
    validate_unitaryExponentIsReal(exponent, __func__);
 
    // notice exponent is a 'qcomp' not a 'qreal' despite validation; see applyMultiStateControlledDiagMatrPower()
 
    // harmlessly re-validates
    applyMultiStateControlledDiagMatrPower(qureg, &control, nullptr, 1, targets, numTargets, matrix, exponent);
}
 
void applyMultiControlledDiagMatrPower(Qureg qureg, int* controls, int numControls, int* targets, int numTargets, DiagMatr matrix, qcomp exponent) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
    validate_matrixDimMatchesTargets(matrix, numTargets, __func__); // also checks fields and is-synced
    validate_matrixIsUnitary(matrix, __func__);
    validate_matrixExpIsNonDiverging(matrix, exponent, __func__);
    validate_unitaryExponentIsReal(exponent, __func__);
 
    // notice exponent is a 'qcomp' not a 'qreal' despite validation; see applyMultiStateControlledDiagMatrPower()
 
    // harmlessly re-validates
    applyMultiStateControlledDiagMatrPower(qureg, controls, nullptr, numControls, targets, numTargets, matrix, exponent);
}
 
void applyMultiStateControlledDiagMatrPower(Qureg qureg, int* controls, int* states, int numControls, int* targets, int numTargets, DiagMatr matrix, qcomp exponent) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
    validate_controlStates(states, numControls, __func__); // can be nullptr, ignoring numControls
    validate_matrixDimMatchesTargets(matrix, numTargets, __func__);
    validate_matrixIsUnitary(matrix, __func__); 
    validate_matrixExpIsNonDiverging(matrix, exponent, __func__);
    validate_unitaryExponentIsReal(exponent, __func__);
 
    // unlike calcExpecDiagMatrPower which accepts a 'qreal' exponent,
    // this function accepts a 'qcomp'. This is because it always
    // evaluates the relatively numerically unstable pow(qcomp,qcomp)
    // overload, rather than pow(qreal,qreal), since there is no reason
    // to think matrix is real (unitarity permits it to be complex).
    // As such, despite unitarity requiring exponent is also real, we
    // accept a qcomp types so that complex exponents can be passed
    // when numerical validation is disabled without a separate func.
 
    bool conj = false;
    auto ctrlList = lists_getList64(controls, numControls);
    auto stateList = util_getList64OrAllOnes(states, numControls);
    auto targList = lists_getList64(targets,  numTargets);
    localiser_statevec_anyCtrlAnyTargDiagMatr(qureg, ctrlList, stateList, targList, matrix, exponent, conj);
 
    if (!qureg.isDensityMatrix)
        return;
 
    conj = true;
    ctrlList = util_getBraQubits(ctrlList, qureg);
    targList = util_getBraQubits(targList, qureg);
    localiser_statevec_anyCtrlAnyTargDiagMatr(qureg, ctrlList, stateList, targList, matrix, exponent, conj);
}
 
} // end de-mangler
 
void applyDiagMatrPower(Qureg qureg, vector<int> targets, DiagMatr matrix, qcomp exponent) {
 
    applyDiagMatrPower(qureg, targets.data(), targets.size(), matrix, exponent);
}
 
void applyControlledDiagMatrPower(Qureg qureg, int control, vector<int> targets, DiagMatr matrix, qcomp exponent) {
 
    applyControlledDiagMatrPower(qureg, control, targets.data(), targets.size(), matrix, exponent);
}
 
void applyMultiControlledDiagMatrPower(Qureg qureg, vector<int> controls, vector<int> targets, DiagMatr matrix, qcomp exponent) {
 
    applyMultiControlledDiagMatrPower(qureg, controls.data(), controls.size(), targets.data(), targets.size(), matrix, exponent);
}
 
void applyMultiStateControlledDiagMatrPower(Qureg qureg, vector<int> controls, vector<int> states, vector<int> targets, DiagMatr matrix, qcomp exponent) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledDiagMatrPower(qureg, controls.data(), states.data(), controls.size(), targets.data(), targets.size(), matrix, exponent);
}
 
 
 
/*
 * FullStateDiagMatr (and power)
 */
 
extern "C" {
 
void applyFullStateDiagMatr(Qureg qureg, FullStateDiagMatr matrix) {
    validate_quregFields(qureg, __func__);
    validate_matrixFields(matrix, __func__);
    validate_matrixAndQuregAreCompatible(matrix, qureg, false, __func__);
    validate_matrixIsUnitary(matrix, __func__);
 
    applyFullStateDiagMatrPower(qureg, matrix, 1); // harmlessly re-validates
}
 
void applyFullStateDiagMatrPower(Qureg qureg, FullStateDiagMatr matrix, qcomp exponent) {
    validate_quregFields(qureg, __func__);
    validate_matrixFields(matrix, __func__);
    validate_matrixAndQuregAreCompatible(matrix, qureg, false, __func__);
    validate_matrixIsUnitary(matrix, __func__);
    validate_unitaryExponentIsReal(exponent, __func__);
    validate_matrixExpIsNonDiverging(matrix, exponent, __func__);
 
    // rho -> matrix^exponent rho conj(matrix^exponent)
    bool leftMultiply = true;
    bool rightMultiply = true;
    bool rightConj = true;
 
    (qureg.isDensityMatrix)?
        localiser_densmatr_allTargDiagMatr(qureg, matrix, exponent, leftMultiply, rightMultiply, rightConj):
        localiser_statevec_allTargDiagMatr(qureg, matrix, exponent);
}
 
} // end de-mangler
 
 
 
/*
 * S gate
 */
 
extern "C" {
 
void applyS(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledS(qureg, nullptr, nullptr, 0, target);
}
 
void applyControlledS(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledS(qureg, &control, nullptr, 1, target);
}
 
void applyMultiControlledS(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledS(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiStateControlledS(Qureg qureg, int* controls, int* states, int numControls, int target) {
 
    static const DiagMatr1 matr = getDiagMatr1({1, 1_i});
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matr, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledS(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledS(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiStateControlledS(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledS(qureg, controls.data(), states.data(), controls.size(), target);
}
 
 
 
/*
 * T gate
 */
 
extern "C" {
 
void applyT(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledT(qureg, nullptr, nullptr, 0, target);
}
 
void applyControlledT(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledT(qureg, &control, nullptr, 1, target);
}
 
void applyMultiControlledT(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledT(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiStateControlledT(Qureg qureg, int* controls, int* states, int numControls, int target) {
 
    static const DiagMatr1 matr = getDiagMatr1({1, (1 + 1_i)/std::sqrt(2)});
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matr, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledT(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledT(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiStateControlledT(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledT(qureg, controls.data(), states.data(), controls.size(), target);
}
 
 
 
/*
 * Hadamard 
 */
 
extern "C" {
 
void applyHadamard(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledHadamard(qureg, nullptr, nullptr, 0, target);
}
 
void applyControlledHadamard(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledHadamard(qureg, &control, nullptr, 1, target);
}
 
void applyMultiControlledHadamard(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledHadamard(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiStateControlledHadamard(Qureg qureg, int* controls, int* states, int numControls, int target) {
 
    static const qcomp a = 1 / std::sqrt(2);
    static const CompMatr1 matr = getCompMatr1({
        {a,  a}, 
        {a, -a}
    });
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matr, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledHadamard(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledHadamard(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiStateControlledHadamard(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledHadamard(qureg, controls.data(), states.data(), controls.size(), target);
}
 
 
 
/*
 * swap
 */
 
extern "C" {
 
void applySwap(Qureg qureg, int qubit1, int qubit2) {
    validate_quregFields(qureg, __func__);
    validate_twoTargets(qureg, qubit1, qubit2, __func__);
 
    // harmlessly re-valdiates
    applyMultiStateControlledSwap(qureg, nullptr, nullptr, 0, qubit1, qubit2);
}
 
void applyControlledSwap(Qureg qureg, int control, int qubit1, int qubit2) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTwoTargets(qureg, control, qubit1, qubit2, __func__);
 
    // harmlessly re-valdiates
    applyMultiStateControlledSwap(qureg, &control, nullptr, 1, qubit1, qubit2);
}
 
void applyMultiControlledSwap(Qureg qureg, int* controls, int numControls, int qubit1, int qubit2) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTwoTargets(qureg, controls, numControls, qubit1, qubit2, __func__);
 
    // harmlessly re-valdiates
    applyMultiStateControlledSwap(qureg, controls, nullptr, numControls, qubit1, qubit2);
}
 
void applyMultiStateControlledSwap(Qureg qureg, int* controls, int* states, int numControls, int qubit1, int qubit2) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTwoTargets(qureg, controls, numControls, qubit1, qubit2, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    auto ctrlList = lists_getList64(controls, numControls);
    auto stateList = util_getList64OrAllOnes(states, numControls);
    localiser_statevec_anyCtrlSwap(qureg, ctrlList, stateList, qubit1, qubit2);
 
    if (!qureg.isDensityMatrix)
        return;
 
    ctrlList = util_getBraQubits(ctrlList, qureg);
    qubit1 = util_getBraQubit(qubit1, qureg);
    qubit2 = util_getBraQubit(qubit2, qureg);
    localiser_statevec_anyCtrlSwap(qureg, ctrlList, stateList, qubit1, qubit2);
}
 
} // end de-mangler
 
void applyMultiControlledSwap(Qureg qureg, vector<int> controls, int qubit1, int qubit2) {
 
    applyMultiControlledSwap(qureg, controls.data(), controls.size(), qubit1, qubit2);
}
 
void applyMultiStateControlledSwap(Qureg qureg, vector<int> controls, vector<int> states, int qubit1, int qubit2) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledSwap(qureg, controls.data(), states.data(), controls.size(), qubit1, qubit2);
}
 
 
 
/*
 * sqrt swap
 */
 
extern "C" {
 
void applySqrtSwap(Qureg qureg, int target1, int target2) {
    validate_quregFields(qureg, __func__);
    validate_twoTargets(qureg,target1, target2, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledSqrtSwap(qureg, nullptr, nullptr, 0, target1, target2);
}
 
void applyControlledSqrtSwap(Qureg qureg, int control, int target1, int target2) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTwoTargets(qureg, control, target1, target2, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledSqrtSwap(qureg, &control, nullptr, 1, target1, target2);
}
 
void applyMultiControlledSqrtSwap(Qureg qureg, int* controls, int numControls, int target1, int target2) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTwoTargets(qureg, controls, numControls, target1, target2, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledSqrtSwap(qureg, controls, nullptr, numControls, target1, target2);
}
 
void applyMultiStateControlledSqrtSwap(Qureg qureg, int* controls, int* states, int numControls, int target1, int target2) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTwoTargets(qureg, controls, numControls, target1, target2, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    /// @todo
    /// this function exacts sqrtSwap as a dense 2-qubit matrix,
    /// where as bespoke communication and simulation strategy is
    /// clearly possible which we have not supported because the gate
    /// is somewhat esoteric. As such, we must validate mixed-amps fit
 
    validate_mixedAmpsFitInNode(qureg, 2, __func__); // to throw SqrtSwap error, not generic CompMatr2 error
 
    static const CompMatr2 matr = getCompMatr2({
        {1, 0, 0, 0},
        {0, .5+.5_i, .5-.5_i, 0},
        {0, .5-.5_i, .5+.5_i, 0},
        {0, 0, 0, 1}});
 
    applyMultiStateControlledCompMatr2(qureg, controls, states, numControls, target1, target2, matr);
}
 
} // end de-mangler
 
void applyMultiControlledSqrtSwap(Qureg qureg, vector<int> controls, int qubit1, int qubit2) {
 
    applyMultiControlledSqrtSwap(qureg, controls.data(), controls.size(), qubit1, qubit2);
}
 
void applyMultiStateControlledSqrtSwap(Qureg qureg, vector<int> controls, vector<int> states, int qubit1, int qubit2) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledSqrtSwap(qureg, controls.data(), states.data(), controls.size(), qubit1, qubit2);
}
 
 
 
/*
 * individual Paulis
 *
 * where Y and Z are most efficiently effected as DiagMatr1,
 * but where X is best effected as a 1-qubit PauliStr.
 */
 
extern "C" {
 
void applyPauliX(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliX(qureg, nullptr, nullptr, 0, target);
}
 
void applyPauliY(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliY(qureg, nullptr, nullptr, 0, target);
}
 
void applyPauliZ(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliZ(qureg, nullptr, nullptr, 0, target);
}
 
void applyControlledPauliX(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliX(qureg, &control, nullptr, 1, target);
}
 
void applyControlledPauliY(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliY(qureg, &control, nullptr, 1, target);
}
 
void applyControlledPauliZ(Qureg qureg, int control, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliZ(qureg, &control, nullptr, 1, target);
}
 
void applyMultiControlledPauliX(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliX(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiControlledPauliY(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliY(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiControlledPauliZ(Qureg qureg, int* controls, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliZ(qureg, controls, nullptr, numControls, target);
}
 
void applyMultiStateControlledPauliX(Qureg qureg, int* controls, int* states, int numControls, int target)  {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // note that for the single-target scenario, we do not call the backend of
    // applyMultiStateControlledPauliStr() since it contains sub-optimal logic
    // which sees the factor of every amplitude dynamically evaluated (based on
    // index parity, etc); the dense-matrix element lookup is faster
 
    /// @todo
    /// a bespoke all-pauli-X function (like in QuEST v3) will be faster still 
    /// since it avoids all superfluous flops; check worthwhile for multi-qubit
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    static const CompMatr1 matrix = util_getPauliX();
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledPauliY(Qureg qureg, int* controls, int* states, int numControls, int target) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    static const CompMatr1 matrix = util_getPauliY();
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledPauliZ(Qureg qureg, int* controls, int* states, int numControls, int target)  {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    static const DiagMatr1 matrix = util_getPauliZ();
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledPauliX(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledPauliX(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiControlledPauliY(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledPauliY(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiControlledPauliZ(Qureg qureg, vector<int> controls, int target) {
 
    applyMultiControlledPauliZ(qureg, controls.data(), controls.size(), target);
}
 
void applyMultiStateControlledPauliX(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPauliX(qureg, controls.data(), states.data(), controls.size(), target);
}
 
void applyMultiStateControlledPauliY(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPauliY(qureg, controls.data(), states.data(), controls.size(), target);
}
 
void applyMultiStateControlledPauliZ(Qureg qureg, vector<int> controls, vector<int> states, int target) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPauliZ(qureg, controls.data(), states.data(), controls.size(), target);
}
 
 
 
/*
 * Pauli strings
 */
 
extern "C" {
 
void applyPauliStr(Qureg qureg, PauliStr str) {
    validate_quregFields(qureg, __func__);
    validate_pauliStrTargets(qureg, str, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliStr(qureg, nullptr, nullptr, 0, str);
}
 
void applyControlledPauliStr(Qureg qureg, int control, PauliStr str) {
    validate_quregFields(qureg, __func__);
    validate_controlAndPauliStrTargets(qureg, control, str, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliStr(qureg, &control, nullptr, 1, str);
}
 
void applyMultiControlledPauliStr(Qureg qureg, int* controls, int numControls, PauliStr str) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndPauliStrTargets(qureg, controls, numControls, str, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliStr(qureg, controls, nullptr, numControls, str);
}
 
void applyMultiStateControlledPauliStr(Qureg qureg, int* controls, int* states, int numControls, PauliStr str) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndPauliStrTargets(qureg, controls, numControls, str, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    qcomp factor = 1;
    auto ctrlList = lists_getList64(controls, numControls);
    auto stateList = util_getList64OrAllOnes(states, numControls);
 
    // when there are no control qubits, we can merge the density matrix's 
    // operation sinto a single tensor, i.e. +- (shift(str) (x) str), to 
    // avoid superfluous re-enumeration of the state
    if (qureg.isDensityMatrix && numControls == 0) {
        factor = paulis_getSignOfPauliStrConj(str);
        ctrlList = util_getConcatenated(ctrlList, util_getBraQubits(ctrlList, qureg));
        stateList = util_getConcatenated(stateList, stateList); 
        str = paulis_getKetAndBraPauliStr(str, qureg);
    }
 
    localiser_statevec_anyCtrlPauliTensor(qureg, ctrlList, stateList, str, factor);
 
    // but density-matrix control qubits require two distinct operations
    if (qureg.isDensityMatrix && numControls > 0) {
        factor = paulis_getSignOfPauliStrConj(str);
        ctrlList = util_getBraQubits(ctrlList, qureg);
        str = paulis_getShiftedPauliStr(str, qureg.numQubits);
        localiser_statevec_anyCtrlPauliTensor(qureg, ctrlList, stateList, str, factor);
    }
}
 
} // end de-mangler
 
void applyMultiControlledPauliStr(Qureg qureg, vector<int> controls, PauliStr str) {
 
    applyMultiControlledPauliStr(qureg, controls.data(), controls.size(), str);
}
 
void applyMultiStateControlledPauliStr(Qureg qureg, vector<int> controls, vector<int> states, PauliStr str) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPauliStr(qureg, controls.data(), states.data(), controls.size(), str);
}
 
 
 
/*
 * individual axis rotations
 */
 
extern "C" {
 
void applyRotateX(Qureg qureg, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateX(qureg, nullptr, nullptr, 0, target, angle);
}
 
void applyRotateY(Qureg qureg, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateY(qureg, nullptr, nullptr, 0, target, angle);
}
 
void applyRotateZ(Qureg qureg, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateZ(qureg, nullptr, nullptr, 0, target, angle);
}
 
void applyControlledRotateX(Qureg qureg, int control, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateX(qureg, &control, nullptr, 1, target, angle);
}
 
void applyControlledRotateY(Qureg qureg, int control, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateY(qureg, &control, nullptr, 1, target, angle);
}
 
void applyControlledRotateZ(Qureg qureg, int control, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, control, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateZ(qureg, &control, nullptr, 1, target, angle);
}
 
void applyMultiControlledRotateX(Qureg qureg, int* controls, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateX(qureg, controls, nullptr, numControls, target, angle);
}
 
void applyMultiControlledRotateY(Qureg qureg, int* controls, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateY(qureg, controls, nullptr, numControls, target, angle);
}
 
void applyMultiControlledRotateZ(Qureg qureg, int* controls, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledRotateZ(qureg, controls, nullptr, numControls, target, angle);
}
 
void applyMultiStateControlledRotateX(Qureg qureg, int* controls, int* states, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // note that for the single-target scenario, we do not call the backend of
    // applyMultiStateControlledPauliGadget() since it contains sub-optimal logic
    // which sees the factor of every amplitude dynamically evaluated (based on
    // index parity, etc); the dense-matrix element lookup is faster
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    CompMatr1 matrix = util_getExpPauliX(angle);
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledRotateY(Qureg qureg, int* controls, int* states, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // note that for the single-target scenario, we do not call the backend of
    // applyMultiStateControlledPauliGadget() since it contains sub-optimal logic
    // which sees the factor of every amplitude dynamically evaluated (based on
    // index parity, etc); the dense-matrix element lookup is faster
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    CompMatr1 matrix = util_getExpPauliY(angle);
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
void applyMultiStateControlledRotateZ(Qureg qureg, int* controls, int* states, int numControls, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, controls, numControls, target, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // note that for the single-target scenario, we do not call the backend of
    // applyMultiStateControlledPauliGadget() since it contains sub-optimal logic
    // which sees the factor of every amplitude dynamically evaluated (based on
    // index parity, etc); the dense-matrix element lookup is faster
 
    // harmlessly re-validates, including hardcoded matrix unitarity
    DiagMatr1 matrix = util_getExpPauliZ(angle);
    validateAndApplyAnyCtrlAnyTargUnitaryMatrix(qureg, controls, states, numControls, &target, 1, matrix, __func__);
}
 
} // end de-mangler
 
void applyMultiControlledRotateX(Qureg qureg, vector<int> controls, int target, qreal angle) {
 
    applyMultiControlledRotateX(qureg, controls.data(), controls.size(), target, angle);
}
 
void applyMultiControlledRotateY(Qureg qureg, vector<int> controls, int target, qreal angle) {
 
    applyMultiControlledRotateY(qureg, controls.data(), controls.size(), target, angle);
}
 
void applyMultiControlledRotateZ(Qureg qureg, vector<int> controls, int target, qreal angle) {
 
    applyMultiControlledRotateZ(qureg, controls.data(), controls.size(), target, angle);
}
 
void applyMultiStateControlledRotateX(Qureg qureg, vector<int> controls, vector<int> states, int target, qreal angle) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledRotateX(qureg, controls.data(), states.data(), controls.size(), target, angle);
}
 
void applyMultiStateControlledRotateY(Qureg qureg, vector<int> controls, vector<int> states, int target, qreal angle) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledRotateY(qureg, controls.data(), states.data(), controls.size(), target, angle);
}
 
void applyMultiStateControlledRotateZ(Qureg qureg, vector<int> controls, vector<int> states, int target, qreal angle) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledRotateZ(qureg, controls.data(), states.data(), controls.size(), target, angle);
}
 
 
 
/*
 * arbitrary axis rotation
 */
 
extern "C" {
 
void applyRotateAroundAxis(Qureg qureg, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, targ, __func__);
    validate_rotationAxisNotZeroVector(axisX, axisY, axisZ, __func__);
 
    applyMultiStateControlledRotateAroundAxis(qureg, nullptr, nullptr, 0, targ, angle, axisX, axisY, axisZ);
}
 
void applyControlledRotateAroundAxis(Qureg qureg, int ctrl, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTarget(qureg, ctrl, targ, __func__);
    validate_rotationAxisNotZeroVector(axisX, axisY, axisZ, __func__);
 
    applyMultiStateControlledRotateAroundAxis(qureg, &ctrl, nullptr, 1, targ, angle, axisX, axisY, axisZ);
}
 
void applyMultiControlledRotateAroundAxis(Qureg qureg, int* ctrls, int numCtrls, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, ctrls, numCtrls, targ, __func__);
    validate_rotationAxisNotZeroVector(axisX, axisY, axisZ, __func__);
 
    applyMultiStateControlledRotateAroundAxis(qureg, ctrls, nullptr, numCtrls, targ, angle, axisX, axisY, axisZ);
}
 
void applyMultiStateControlledRotateAroundAxis(Qureg qureg, int* ctrls, int* states, int numCtrls, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTarget(qureg, ctrls, numCtrls, targ, __func__);
    validate_controlStates(states, numCtrls, __func__); // permits states==nullptr
    validate_rotationAxisNotZeroVector(axisX, axisY, axisZ, __func__);
 
    // defer division of vector norm to improve numerical accuracy
    qreal norm = std::sqrt(std::pow(axisX,2) + std::pow(axisY,2) + std::pow(axisZ,2)); // != 0
 
    // treat as generic 1-qubit matrix
    qreal c = std::cos(angle/2);
    qreal s = std::sin(angle/2);
    qcomp u11 = c - (s * axisZ * 1_i) / norm;
    qcomp u12 =   - (s * (axisY + axisX * 1_i)) / norm;
    qcomp u21 =     (s * (axisY - axisX * 1_i)) / norm;
    qcomp u22 = c + (s * axisZ * 1_i) / norm;
    auto matr = getCompMatr1({{u11,u12},{u21,u22}});
 
    // harmlessly re-validates, and checks unitarity of matr
    applyMultiStateControlledCompMatr1(qureg, ctrls, states, numCtrls, targ, matr);
}
 
} // end de-mangler
 
void applyMultiControlledRotateAroundAxis(Qureg qureg, vector<int> ctrls, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
 
    applyMultiControlledRotateAroundAxis(qureg, ctrls.data(), ctrls.size(), targ, angle, axisX, axisY, axisZ);
}
 
void applyMultiStateControlledRotateAroundAxis(Qureg qureg, vector<int> ctrls, vector<int> states, int targ, qreal angle, qreal axisX, qreal axisY, qreal axisZ) {
    validate_controlsMatchStates(ctrls.size(), states.size(), __func__);
 
    applyMultiStateControlledRotateAroundAxis(qureg, ctrls.data(), states.data(), ctrls.size(), targ, angle, axisX, axisY, axisZ);
}
 
 
 
/*
 * Pauli gadgets
 */
 
extern "C" {
 
void applyPauliGadget(Qureg qureg, PauliStr str, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_pauliStrTargets(qureg, str, __func__);
    
    applyMultiStateControlledPauliGadget(qureg, nullptr, nullptr, 0, str, angle);
}
 
void applyNonUnitaryPauliGadget(Qureg qureg, PauliStr str, qcomp angle) {
    validate_quregFields(qureg, __func__);
    validate_pauliStrTargets(qureg, str, __func__);
 
    qcomp phase = util_getPhaseFromGateAngle(angle);
    auto none = lists_getEmptyList64();
    localiser_statevec_anyCtrlPauliGadget(qureg, none, none, str, phase);
 
    if (!qureg.isDensityMatrix)
        return;
 
    // conj(e^i(a)P) = e^(-i s conj(a) P)
    phase = - std::conj(phase) * paulis_getSignOfPauliStrConj(str);
    str = paulis_getShiftedPauliStr(str, qureg.numQubits);
    localiser_statevec_anyCtrlPauliGadget(qureg, none, none, str, phase);
}
 
void applyControlledPauliGadget(Qureg qureg, int control, PauliStr str, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlAndPauliStrTargets(qureg, control, str, __func__);
    
    applyMultiStateControlledPauliGadget(qureg, &control, nullptr, 1, str, angle);
}
 
void applyMultiControlledPauliGadget(Qureg qureg, int* controls, int numControls, PauliStr str, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndPauliStrTargets(qureg, controls, numControls, str, __func__);
    
    applyMultiStateControlledPauliGadget(qureg, controls, nullptr, numControls, str, angle);
}
 
void applyMultiStateControlledPauliGadget(Qureg qureg, int* controls, int* states, int numControls, PauliStr str, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndPauliStrTargets(qureg, controls, numControls, str, __func__);
    validate_controlStates(states, numControls, __func__); // permits states==nullptr
 
    // a non-controlled str=I effects a global phase change (of -angle/2) which does not 
    // at all change a density matrix; the subsequent dagger operation would undo it,
    // which we avoid to preserve numerical accuracy
    if (paulis_isIdentity(str) && numControls == 0 && qureg.isDensityMatrix)
        return;
 
    // when numControls >= 1, all amps satisfying the control condition undergo a phase 
    // change of -angle/2, as if all non-control-qubits were targeted by exp(-angle/2)I,
    // which is sufficiently efficient using the existing gadget backend function
 
    qreal phase = util_getPhaseFromGateAngle(angle);
    auto ctrlList = lists_getList64(controls, numControls);
    auto stateList = util_getList64OrAllOnes(states, numControls);
    localiser_statevec_anyCtrlPauliGadget(qureg, ctrlList, stateList, str, phase);
 
    if (!qureg.isDensityMatrix)
        return;
 
    // conj(e^(i a P)) = e^(-i s a P)
    phase *= - paulis_getSignOfPauliStrConj(str);
    ctrlList = util_getBraQubits(ctrlList, qureg);
    str = paulis_getShiftedPauliStr(str, qureg.numQubits);
    localiser_statevec_anyCtrlPauliGadget(qureg, ctrlList, stateList, str, phase);
}
 
} // end de-mangler
 
void applyMultiControlledPauliGadget(Qureg qureg, vector<int> controls, PauliStr str, qreal angle) {
 
    applyMultiControlledPauliGadget(qureg, controls.data(), controls.size(), str, angle);
}
 
void applyMultiStateControlledPauliGadget(Qureg qureg, vector<int> controls, vector<int> states, PauliStr str, qreal angle) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPauliGadget(qureg, controls.data(), states.data(), controls.size(), str, angle);
}
 
 
 
/*
 * phase gadgets
 */
 
extern "C" {
 
void applyPhaseGadget(Qureg qureg, int* targets, int numTargets, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPhaseGadget(qureg, nullptr, nullptr, 0, targets, numTargets, angle);
}
 
void applyControlledPhaseGadget(Qureg qureg, int control, int* targets, int numTargets, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTargets(qureg, control, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPhaseGadget(qureg, &control, nullptr, 1, targets, numTargets, angle);
}
 
void applyMultiControlledPhaseGadget(Qureg qureg, int* controls, int numControls, int* targets, int numTargets, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledPhaseGadget(qureg, controls, nullptr, numControls, targets, numTargets, angle);
}
 
void applyMultiStateControlledPhaseGadget(Qureg qureg, int* controls, int* states, int numControls, int* targets, int numTargets, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
    validate_controlStates(states, numControls, __func__);
 
    qreal phase = util_getPhaseFromGateAngle(angle);
    auto ctrlList = lists_getList64(controls, numControls);
    auto stateList = util_getList64OrAllOnes(states, numControls);
    auto targList = lists_getList64(targets,  numTargets);
    localiser_statevec_anyCtrlPhaseGadget(qureg, ctrlList, stateList, targList, phase);
 
    if (!qureg.isDensityMatrix)
        return;
 
    phase *= -1;
    ctrlList = util_getBraQubits(ctrlList, qureg);
    targList = util_getBraQubits(targList, qureg);
    localiser_statevec_anyCtrlPhaseGadget(qureg, ctrlList, stateList, targList, phase);
}
 
} // end de-mangler
 
void applyPhaseGadget(Qureg qureg, vector<int> targets, qreal angle) {
 
    applyPhaseGadget(qureg, targets.data(), targets.size(), angle);
}
 
void applyControlledPhaseGadget(Qureg qureg, int control, vector<int> targets, qreal angle) {
 
    applyControlledPhaseGadget(qureg, control, targets.data(), targets.size(), angle);
}
 
void applyMultiControlledPhaseGadget(Qureg qureg, vector<int> controls, vector<int> targets, qreal angle) {
 
    applyMultiControlledPhaseGadget(qureg, controls.data(), controls.size(), targets.data(), targets.size(), angle);
}
 
void applyMultiStateControlledPhaseGadget(Qureg qureg, vector<int> controls, vector<int> states, vector<int> targets, qreal angle) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledPhaseGadget(qureg, controls.data(), states.data(), controls.size(), targets.data(), targets.size(), angle);
}
 
 
 
/*
 * phase shift
 */
 
extern "C" {
 
void applyPhaseShift(Qureg qureg, int target, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiQubitPhaseShift(qureg, &target, 1, angle);
}
 
void applyTwoQubitPhaseShift(Qureg qureg, int target1, int target2, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_twoTargets(qureg, target1, target2, __func__);
 
    // harmlessly re-validates
    int targets[] = {target1, target2};
    applyMultiQubitPhaseShift(qureg, targets, 2, angle);
}
 
void applyMultiQubitPhaseShift(Qureg qureg, int* targets, int numTargets, qreal angle) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
 
    // treat as a (numTargets-1)-controlled 1-target diagonal matrix
    static DiagMatr1 matr = getDiagMatr1({1, /*un-init*/ 0});
    matr.elems[1] = std::exp(1_i * angle); // micro-optimisation
 
    // harmlessly re-validates
    applyMultiStateControlledDiagMatr1(qureg, &targets[1], nullptr, numTargets-1, targets[0], matr);
}
 
} // end de-mangler
 
void applyMultiQubitPhaseShift(Qureg qureg, vector<int> targets, qreal angle) {
 
    applyMultiQubitPhaseShift(qureg, targets.data(), targets.size(), angle);
}
 
 
 
/*
 * phase flips
 */
 
extern "C" {
 
void applyPhaseFlip(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // harmlessly re-validates
    applyMultiQubitPhaseFlip(qureg, &target, 1);
}
 
void applyTwoQubitPhaseFlip(Qureg qureg, int target1, int target2) {
    validate_quregFields(qureg, __func__);
    validate_twoTargets(qureg, target1, target2, __func__);
 
    // harmlessly re-validates
    int targets[] = {target1, target2};
    applyMultiQubitPhaseFlip(qureg, targets, 2);
}
 
void applyMultiQubitPhaseFlip(Qureg qureg, int* targets, int numTargets) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
 
    // treat as a (numTargets-1)-controlled 1-target Pauli Z
    static const DiagMatr1 matr = getDiagMatr1({1, -1});
 
    // harmlessly re-validates
    applyMultiStateControlledDiagMatr1(qureg, &targets[1], nullptr, numTargets-1, targets[0], matr);
}
 
} // end de-mangler
 
void applyMultiQubitPhaseFlip(Qureg qureg, vector<int> targets) {
 
    applyMultiQubitPhaseFlip(qureg, targets.data(), targets.size());
}
 
 
 
/*
 * many-qubit NOTs and CNOTs
 */
 
extern "C" {
 
void applyMultiQubitNot(Qureg qureg, int* targets, int numTargets) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledMultiQubitNot(qureg, nullptr, nullptr, 0, targets, numTargets);
}
 
void applyControlledMultiQubitNot(Qureg qureg, int control, int* targets, int numTargets) {
    validate_quregFields(qureg, __func__);
    validate_controlAndTargets(qureg, control, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledMultiQubitNot(qureg, &control, nullptr, 1, targets, numTargets);
}
 
void applyMultiControlledMultiQubitNot(Qureg qureg, int* controls, int numControls, int* targets, int numTargets) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
 
    // harmlessly re-validates
    applyMultiStateControlledMultiQubitNot(qureg, controls, nullptr, numControls, targets, numTargets);
}
 
void applyMultiStateControlledMultiQubitNot(Qureg qureg, int* controls, int* states, int numControls, int* targets, int numTargets) {
    validate_quregFields(qureg, __func__);
    validate_controlsAndTargets(qureg, controls, numControls, targets, numTargets, __func__);
    validate_controlStates(states, numControls, __func__);
 
    // treat as an all-X PauliStr
    PauliStr str = getPauliStr(std::string(numTargets, 'X'), targets, numTargets);
 
    // harmlessly re-validates
    applyMultiStateControlledPauliStr(qureg, controls, states, numControls, str);
}
 
} // end de-mangler
 
void applyMultiQubitNot(Qureg qureg, vector<int> targets) {
 
    applyMultiQubitNot(qureg, targets.data(), targets.size());
}
 
void applyControlledMultiQubitNot(Qureg qureg, int control, vector<int> targets) {
 
    applyControlledMultiQubitNot(qureg, control, targets.data(), targets.size());
}
 
void applyMultiControlledMultiQubitNot(Qureg qureg, vector<int> controls, vector<int> targets) {
 
    applyMultiControlledMultiQubitNot(qureg, controls.data(), controls.size(), targets.data(), targets.size());
}
 
void applyMultiStateControlledMultiQubitNot(Qureg qureg, vector<int> controls, vector<int> states, vector<int> targets) {
    validate_controlsMatchStates(controls.size(), states.size(), __func__);
 
    applyMultiStateControlledMultiQubitNot(qureg, controls.data(), states.data(), controls.size(), targets.data(), targets.size());
}
 
 
 
/*
 * projectors
 */
 
extern "C" {
 
void applyQubitProjector(Qureg qureg, int target, int outcome) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
    validate_measurementOutcomeIsValid(outcome, __func__); 
    
    qreal prob = 1;
 
    auto targList    = lists_getList64({target});
    auto outcomeList = lists_getList64({outcome});
 
    // density matrix has an optimised func in lieu of calling the statevector func twice
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, targList, outcomeList, prob):
        localiser_statevec_multiQubitProjector(qureg, targList, outcomeList, prob);
}
 
void applyMultiQubitProjector(Qureg qureg, int* qubits, int* outcomes, int numQubits) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, qubits, numQubits, __func__);
    validate_measurementOutcomesAreValid(outcomes, numQubits, __func__);
 
    qreal prob = 1;
    auto qubitList = lists_getList64(qubits, numQubits);
    auto outcomeList = lists_getList64(outcomes, numQubits);
 
    // density matrix has an optimised func in lieu of calling the statevector func twice
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, qubitList, outcomeList, prob):
        localiser_statevec_multiQubitProjector(qureg, qubitList, outcomeList, prob);
}
 
} // end de-mangler
 
void applyMultiQubitProjector(Qureg qureg, vector<int> qubits, vector<int> outcomes) {
    validate_measurementOutcomesMatchTargets(qubits.size(), outcomes.size(), __func__);
 
    applyMultiQubitProjector(qureg, qubits.data(), outcomes.data(), outcomes.size());
}
 
 
 
/*
 * measurement
 */
 
extern "C" {
 
int applyQubitMeasurement(Qureg qureg, int target) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    qreal prob = 0; // ignored
    return applyQubitMeasurementAndGetProb(qureg, target, &prob); // harmlessly re-validates
}
 
int applyQubitMeasurementAndGetProb(Qureg qureg, int target, qreal* probability) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
 
    // we do not assume state normalisation (that is posteriori checked),
    // so we must perform two reductions; one for each outcome. We choose
    // to re-enumerate the state (potentially doubling caching costs) to
    // avoid the nuisances/race-cons of parallel adding to two scalars.
    vector<qreal> probs(2);
    probs[0] = calcProbOfQubitOutcome(qureg, target, 0); // harmlessly re-validates
    probs[1] = calcProbOfQubitOutcome(qureg, target, 1); // " "
    validate_measurementProbsAreNormalised(probs, __func__);
 
    // randomly choose the outcome
    int outcome = rand_getRandomSingleQubitOutcome(probs[0]);
    *probability = probs[outcome];
 
    auto targList    = lists_getList64({target});
    auto outcomeList = lists_getList64({outcome});
 
    // collapse to the outcome
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, targList, outcomeList, *probability):
        localiser_statevec_multiQubitProjector(qureg, targList, outcomeList, *probability);
 
    return outcome;
}
 
qreal applyForcedQubitMeasurement(Qureg qureg, int target, int outcome) {
    validate_quregFields(qureg, __func__);
    validate_target(qureg, target, __func__);
    validate_measurementOutcomeIsValid(outcome, __func__);
 
    // note that we do not merely invoke applyForcedMultiQubitMeasurement()
    // because we must validate the renormalising probability and
    // report this function's name during the error message
    qreal prob = calcProbOfQubitOutcome(qureg, target, outcome); // harmlessly re-validates
    validate_measurementOutcomeProbNotZero(outcome, prob, __func__);
 
    auto targList    = lists_getList64({target});
    auto outcomeList = lists_getList64({outcome});
 
    // project to the outcome, renormalising the surviving states
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, targList, outcomeList, prob):
        localiser_statevec_multiQubitProjector(qureg, targList, outcomeList, prob);
 
    return prob;
}
 
qindex applyMultiQubitMeasurement(Qureg qureg, int* qubits, int numQubits) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, qubits, numQubits, __func__);
 
    qreal prob = 0; // ignored
 
    // below validates post-measurement and would report 'AndGetProb' function suffix. Eh!
    return applyMultiQubitMeasurementAndGetProb(qureg, qubits, numQubits, &prob);
}
 
qindex applyMultiQubitMeasurementAndGetProb(Qureg qureg, int* qubits, int numQubits, qreal* probability) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, qubits, numQubits, __func__);
 
    // find the probability of all possible outcomes...
    qindex numProbs = powerOf2(numQubits);
 
    // by allocating a temp vector, and validating successful (since exponentially big!)
    vector<qreal> probs;
    auto callback = [&]() { validate_tempListAllocSucceeded(false, numProbs, sizeof(qreal), __func__); };
    util_tryAllocVector(probs, numProbs, callback);
 
    // populate probs
    calcProbsOfAllMultiQubitOutcomes(probs.data(), qureg, qubits, numQubits); // harmlessly re-validates
 
    // we cannot meaningfully sample these probs if not normalised
    validate_measurementProbsAreNormalised(probs, __func__);
 
    // randomly choose an outcome
    qindex outcome = rand_getRandomMultiQubitOutcome(probs);
    *probability = probs[outcome];
 
    // map outcome to individual qubit outcomes
    auto qubitList = lists_getList64(qubits, numQubits);
    auto outcomeList = util_getConstantList(-1, numQubits);
    setToBitsOfInteger(outcomeList.data(), outcome, numQubits);
 
    // project to the outcomes, renormalising the surviving states
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, qubitList, outcomeList, *probability):
        localiser_statevec_multiQubitProjector(qureg, qubitList, outcomeList, *probability);
 
    return outcome;
}
 
qreal applyForcedMultiQubitMeasurement(Qureg qureg, int* qubits, int* outcomes, int numQubits) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, qubits, numQubits, __func__);
    validate_measurementOutcomesAreValid(outcomes, numQubits, __func__);
 
    auto qubitList = lists_getList64(qubits, numQubits);
    auto outcomeList = lists_getList64(outcomes, numQubits);
 
    // ensure probability of the forced measurement outcome is not negligible
    qreal prob = calcProbOfMultiQubitOutcome(qureg, qubits, outcomes, numQubits); // harmlessly re-validates
    validate_measurementOutcomesProbNotZero(outcomes, numQubits, prob, __func__);
 
    // project to the outcome, renormalising the surviving states
    (qureg.isDensityMatrix)?
        localiser_densmatr_multiQubitProjector(qureg, qubitList, outcomeList, prob):
        localiser_statevec_multiQubitProjector(qureg, qubitList, outcomeList, prob);
 
    return prob;
}
 
} // end de-mangler
 
qindex applyMultiQubitMeasurement(Qureg qureg, vector<int> qubits) {
 
    return applyMultiQubitMeasurement(qureg, qubits.data(), qubits.size());
}
 
qindex applyMultiQubitMeasurementAndGetProb(Qureg qureg, vector<int> qubits, qreal* probability) {
 
    return applyMultiQubitMeasurementAndGetProb(qureg, qubits.data(), qubits.size(), probability);
}
 
qreal applyForcedMultiQubitMeasurement(Qureg qureg, vector<int> qubits, vector<int> outcomes) {
    validate_measurementOutcomesMatchTargets(qubits.size(), outcomes.size(), __func__);
 
    return applyForcedMultiQubitMeasurement(qureg, qubits.data(), outcomes.data(), outcomes.size());
}
 
 
 
/*
 * QFT
 */
 
extern "C" {
 
void applyQuantumFourierTransform(Qureg qureg, int* targets, int numTargets, bool inverse) {
    validate_quregFields(qureg, __func__);
    validate_targets(qureg, targets, numTargets, __func__);
 
    /// @todo
    /// change this placeholder implementation to the bespoke, optimised routine,
    /// wherein each contiguous controlled-phase gate is merged
 
    int mid = numTargets/2; // floors
 
    // don't think it is worth overcomplicating this to avoid slight repetition
    if (inverse) {
        for (int n=0; n<mid; n++)
            applySwap(qureg, targets[n], targets[numTargets-1-n]);
 
        for (int n=0; n<numTargets; n++) {
            for (int m=0; m<n; m++) {
                qreal arg = - const_PI / powerOf2(m+1);
                applyTwoQubitPhaseShift(qureg, targets[n], targets[n-m-1], arg);
            }
            applyHadamard(qureg, targets[n]);
        }
 
    } else {
        for (int n=numTargets-1; n>=0; n--) {
            applyHadamard(qureg, targets[n]);
 
            for (int m=0; m<n; m++) {
                qreal arg = const_PI / powerOf2(m+1);
                applyTwoQubitPhaseShift(qureg, targets[n], targets[n-m-1], arg);
            }
        }
 
        for (int n=0; n<mid; n++)
            applySwap(qureg, targets[n], targets[numTargets-1-n]);
    }
}
 
void applyFullQuantumFourierTransform(Qureg qureg, bool inverse) {
    validate_quregFields(qureg, __func__);
 
    auto targets = util_getRange(qureg.numQubits);
    applyQuantumFourierTransform(qureg, targets.data(), targets.size(), inverse);
}
 
} // end de-mangler
 
void applyQuantumFourierTransform(Qureg qureg, vector<int> targets, bool inverse) {
 
    applyQuantumFourierTransform(qureg, targets.data(), targets.size(), inverse);
}