Skip to content

Commit 2a7a7ba

Browse files
bmhowe23bmhowe23
committed
Implement new noise models
Signed-off-by: Ben Howe <[email protected]>
1 parent c47d92c commit 2a7a7ba

File tree

6 files changed

+378
-7
lines changed

6 files changed

+378
-7
lines changed

‎runtime/common/NoiseModel.cpp

+11-2
Original file line numberDiff line numberDiff line change
@@ -379,10 +379,19 @@ noise_model::get_channels(const std::string &quantumOp,
379379
}
380380

381381
noise_model::noise_model() {
382-
register_channel<bit_flip_channel>();
383-
register_channel<phase_flip_channel>();
384382
register_channel<depolarization_channel>();
385383
register_channel<amplitude_damping_channel>();
384+
register_channel<bit_flip_channel>();
385+
register_channel<phase_flip_channel>();
386+
register_channel<x_error>();
387+
register_channel<y_error>();
388+
register_channel<z_error>();
389+
register_channel<amplitude_damping>();
390+
register_channel<phase_damping>();
391+
register_channel<pauli1>();
392+
register_channel<pauli2>();
393+
register_channel<depolarization1>();
394+
register_channel<depolarization2>();
386395
}
387396

388397
std::string get_noise_model_type_name(noise_model_type type) {

‎runtime/common/NoiseModel.h

+261
Original file line numberDiff line numberDiff line change
@@ -19,7 +19,43 @@
1919

2020
namespace cudaq::details {
2121
void warn(const std::string_view msg);
22+
23+
/// @brief Typedef for a matrix wrapper using std::vector<cudaq::complex>
24+
using matrix_wrapper = std::vector<cudaq::complex>;
25+
26+
/// @brief Compute the Kronecker product of two matrices
27+
///
28+
/// @param A First matrix
29+
/// @param rowsA Number of rows in matrix A
30+
/// @param colsA Number of columns in matrix A
31+
/// @param B Second matrix
32+
/// @param rowsB Number of rows in matrix B
33+
/// @param colsB Number of columns in matrix B
34+
/// @return matrix_wrapper Result of the Kronecker product
35+
inline matrix_wrapper kron(const matrix_wrapper &A, int rowsA, int colsA,
36+
const matrix_wrapper &B, int rowsB, int colsB) {
37+
matrix_wrapper C((rowsA * rowsB) * (colsA * colsB));
38+
for (int i = 0; i < rowsA; ++i) {
39+
for (int j = 0; j < colsA; ++j) {
40+
for (int k = 0; k < rowsB; ++k) {
41+
for (int l = 0; l < colsB; ++l) {
42+
C[(i * rowsB + k) * (colsA * colsB) + (j * colsB + l)] =
43+
A[i * colsA + j] * B[k * colsB + l];
44+
}
45+
}
46+
}
47+
}
48+
return C;
49+
}
50+
51+
inline matrix_wrapper scale(const cudaq::real s, const matrix_wrapper &A) {
52+
matrix_wrapper result;
53+
result.reserve(A.size());
54+
for (auto a : A)
55+
result.push_back(s * a);
56+
return result;
2257
}
58+
} // namespace cudaq::details
2359

2460
namespace cudaq {
2561

@@ -597,4 +633,229 @@ class phase_flip_channel : public kraus_channel {
597633
REGISTER_KRAUS_CHANNEL(
598634
noise_model_strings[(int)noise_model_type::phase_flip_channel])
599635
};
636+
637+
/// @brief amplitude_damping is the same as amplitude_damping_channel.
638+
class amplitude_damping : public amplitude_damping_channel {
639+
public:
640+
amplitude_damping(const std::vector<cudaq::real> &p)
641+
: amplitude_damping_channel(p) {
642+
noise_type = noise_model_type::amplitude_damping;
643+
}
644+
amplitude_damping(const real probability)
645+
: amplitude_damping_channel(probability) {
646+
noise_type = noise_model_type::amplitude_damping;
647+
}
648+
REGISTER_KRAUS_CHANNEL(
649+
noise_model_strings[(int)noise_model_type::amplitude_damping])
650+
};
651+
652+
/// @brief phase_damping is the same as phase_flip_channel.
653+
class phase_damping : public phase_flip_channel {
654+
public:
655+
phase_damping(const std::vector<cudaq::real> &p) : phase_flip_channel(p) {
656+
noise_type = noise_model_type::phase_damping;
657+
}
658+
phase_damping(const real probability) : phase_flip_channel(probability) {
659+
noise_type = noise_model_type::phase_damping;
660+
}
661+
REGISTER_KRAUS_CHANNEL(
662+
noise_model_strings[(int)noise_model_type::phase_damping])
663+
};
664+
665+
/// @brief z_error is the same as phase_flip_channel.
666+
class z_error : public phase_flip_channel {
667+
public:
668+
z_error(const std::vector<cudaq::real> &p) : phase_flip_channel(p) {
669+
noise_type = noise_model_type::z_error;
670+
}
671+
z_error(const real probability) : phase_flip_channel(probability) {
672+
noise_type = noise_model_type::z_error;
673+
}
674+
REGISTER_KRAUS_CHANNEL(noise_model_strings[(int)noise_model_type::z_error])
675+
};
676+
677+
/// @brief x_error is the same as bit_flip_channel
678+
class x_error : public bit_flip_channel {
679+
public:
680+
x_error(const std::vector<cudaq::real> &p) : bit_flip_channel(p) {
681+
noise_type = noise_model_type::x_error;
682+
}
683+
x_error(const real probability) : bit_flip_channel(probability) {
684+
noise_type = noise_model_type::x_error;
685+
}
686+
REGISTER_KRAUS_CHANNEL(noise_model_strings[(int)noise_model_type::x_error])
687+
};
688+
689+
/// @brief y_error is the same as bit_flip_channel
690+
class y_error : public kraus_channel {
691+
public:
692+
constexpr static std::size_t num_parameters = 1;
693+
constexpr static std::size_t num_targets = 1;
694+
y_error(const std::vector<cudaq::real> &p) {
695+
cudaq::real probability = p[0];
696+
std::complex<cudaq::real> i{0, 1};
697+
std::vector<cudaq::complex> k0v{std::sqrt(1 - probability), 0, 0,
698+
std::sqrt(1 - probability)},
699+
k1v{0, -i * std::sqrt(probability), i * std::sqrt(probability), 0};
700+
ops = {k0v, k1v};
701+
this->parameters.push_back(probability);
702+
noise_type = noise_model_type::y_error;
703+
validateCompleteness();
704+
generateUnitaryParameters();
705+
}
706+
y_error(const real probability)
707+
: y_error(std::vector<cudaq::real>{probability}) {}
708+
REGISTER_KRAUS_CHANNEL(noise_model_strings[(int)noise_model_type::y_error])
709+
};
710+
711+
class pauli1 : public kraus_channel {
712+
public:
713+
constexpr static std::size_t num_parameters = 3;
714+
constexpr static std::size_t num_targets = 1;
715+
pauli1(const std::vector<cudaq::real> &p) {
716+
if (p.size() == num_parameters)
717+
throw std::runtime_error(
718+
"Invalid number of elements to pauli1 constructor. Must be 3.");
719+
cudaq::real sum = 0;
720+
for (auto pp : p)
721+
sum += pp;
722+
// This is just a first-level error check. Additional checks are done in
723+
// validateCompleteness.
724+
if (sum > static_cast<cudaq::real>(1.0 + 1e-6))
725+
throw std::runtime_error("Sum of pauli1 parameters is >1. Must be <= 1.");
726+
727+
std::complex<cudaq::real> i{0, 1};
728+
cudaq::details::matrix_wrapper I({1, 0, 0, 1});
729+
cudaq::details::matrix_wrapper X({0, 1, 1, 0});
730+
cudaq::details::matrix_wrapper Y({0, -i, i, 0});
731+
cudaq::details::matrix_wrapper Z({1, 0, 0, -1});
732+
cudaq::real p0 =
733+
std::sqrt(std::max(static_cast<cudaq::real>(1.0 - p[0] - p[1] - p[2]),
734+
static_cast<cudaq::real>(0)));
735+
cudaq::real px = std::sqrt(p[0]);
736+
cudaq::real py = std::sqrt(p[1]);
737+
cudaq::real pz = std::sqrt(p[2]);
738+
std::vector<cudaq::complex> k0v = details::scale(p0, I);
739+
std::vector<cudaq::complex> k1v = details::scale(px, X);
740+
std::vector<cudaq::complex> k2v = details::scale(py, Y);
741+
std::vector<cudaq::complex> k3v = details::scale(pz, Z);
742+
ops = {k0v, k1v, k2v, k3v};
743+
this->parameters.reserve(p.size());
744+
for (auto pp : p)
745+
this->parameters.push_back(pp);
746+
noise_type = cudaq::noise_model_type::pauli1;
747+
validateCompleteness();
748+
generateUnitaryParameters();
749+
}
750+
REGISTER_KRAUS_CHANNEL()
751+
};
752+
753+
class pauli2 : public kraus_channel {
754+
public:
755+
constexpr static std::size_t num_parameters = 15;
756+
constexpr static std::size_t num_targets = 2;
757+
pauli2(const std::vector<cudaq::real> &p) {
758+
if (p.size() == num_parameters)
759+
throw std::runtime_error(
760+
"Invalid number of elements to pauli2 constructor. Must be 15.");
761+
cudaq::real sum = 0;
762+
for (auto pp : p)
763+
sum += pp;
764+
// This is just a first-level error check. Additional checks are done in
765+
// validateCompleteness.
766+
if (sum > static_cast<cudaq::real>(1.0 + 1e-6))
767+
throw std::runtime_error("Sum of pauli2 parameters is >1. Must be <= 1.");
768+
769+
std::complex<cudaq::real> i{0, 1};
770+
cudaq::details::matrix_wrapper I({1, 0, 0, 1});
771+
cudaq::details::matrix_wrapper X({0, 1, 1, 0});
772+
cudaq::details::matrix_wrapper Y({0, -i, i, 0});
773+
cudaq::details::matrix_wrapper Z({1, 0, 0, -1});
774+
cudaq::real pii = std::sqrt(std::max(static_cast<cudaq::real>(1.0 - sum),
775+
static_cast<cudaq::real>(0)));
776+
777+
ops.push_back(details::scale(pii, details::kron(I, 2, 2, I, 2, 2)));
778+
ops.push_back(details::scale(p[0], details::kron(I, 2, 2, X, 2, 2)));
779+
ops.push_back(details::scale(p[1], details::kron(I, 2, 2, Y, 2, 2)));
780+
ops.push_back(details::scale(p[2], details::kron(I, 2, 2, Z, 2, 2)));
781+
782+
ops.push_back(details::scale(p[3], details::kron(X, 2, 2, I, 2, 2)));
783+
ops.push_back(details::scale(p[4], details::kron(X, 2, 2, X, 2, 2)));
784+
ops.push_back(details::scale(p[5], details::kron(X, 2, 2, Y, 2, 2)));
785+
ops.push_back(details::scale(p[6], details::kron(X, 2, 2, Z, 2, 2)));
786+
787+
ops.push_back(details::scale(p[7], details::kron(Y, 2, 2, I, 2, 2)));
788+
ops.push_back(details::scale(p[8], details::kron(Y, 2, 2, X, 2, 2)));
789+
ops.push_back(details::scale(p[9], details::kron(Y, 2, 2, Y, 2, 2)));
790+
ops.push_back(details::scale(p[10], details::kron(Y, 2, 2, Z, 2, 2)));
791+
792+
ops.push_back(details::scale(p[11], details::kron(Z, 2, 2, I, 2, 2)));
793+
ops.push_back(details::scale(p[12], details::kron(Z, 2, 2, X, 2, 2)));
794+
ops.push_back(details::scale(p[13], details::kron(Z, 2, 2, Y, 2, 2)));
795+
ops.push_back(details::scale(p[14], details::kron(Z, 2, 2, Z, 2, 2)));
796+
797+
this->parameters.reserve(p.size());
798+
for (auto pp : p)
799+
this->parameters.push_back(pp);
800+
noise_type = cudaq::noise_model_type::pauli2;
801+
validateCompleteness();
802+
generateUnitaryParameters();
803+
}
804+
REGISTER_KRAUS_CHANNEL()
805+
};
806+
807+
/// @brief depolarization1 is the same as depolarization_channel
808+
class depolarization1 : public depolarization_channel {
809+
public:
810+
depolarization1(const std::vector<cudaq::real> &p)
811+
: depolarization_channel(p) {
812+
noise_type = noise_model_type::depolarization1;
813+
}
814+
depolarization1(const real probability)
815+
: depolarization_channel(probability) {
816+
noise_type = noise_model_type::depolarization1;
817+
}
818+
REGISTER_KRAUS_CHANNEL()
819+
};
820+
821+
class depolarization2 : public kraus_channel {
822+
public:
823+
constexpr static std::size_t num_parameters = 1;
824+
constexpr static std::size_t num_targets = 2;
825+
depolarization2(const std::vector<cudaq::real> p) : kraus_channel() {
826+
auto three = static_cast<cudaq::real>(3.);
827+
auto negOne = static_cast<cudaq::real>(-1.);
828+
auto probability = p[0];
829+
830+
std::vector<std::vector<cudaq::complex>> singleQubitKraus = {
831+
{std::sqrt(1 - probability), 0, 0, std::sqrt(1 - probability)},
832+
{0, std::sqrt(probability / three), std::sqrt(probability / three), 0},
833+
{0, cudaq::complex{0, negOne * std::sqrt(probability / three)},
834+
cudaq::complex{0, std::sqrt(probability / three)}, 0},
835+
{std::sqrt(probability / three), 0, 0,
836+
negOne * std::sqrt(probability / three)}};
837+
838+
// Generate 2-qubit Kraus operators
839+
for (const auto &k1 : singleQubitKraus) {
840+
for (const auto &k2 : singleQubitKraus) {
841+
ops.push_back(details::kron(k1, 2, 2, k2, 2, 2));
842+
}
843+
}
844+
this->parameters.push_back(probability);
845+
noise_type = cudaq::noise_model_type::depolarization2;
846+
validateCompleteness();
847+
generateUnitaryParameters();
848+
}
849+
850+
/// @brief Construct a two qubit kraus channel that applies a depolarization
851+
/// channel on either qubit independently.
852+
///
853+
/// @param p The probability of any depolarizing error happening in the 2
854+
/// qubits. (Setting this to 1.0 ensures that "II" cannot happen; maximal
855+
/// mixing occurs at p = 0.9375.)
856+
depolarization2(const real probability)
857+
: depolarization2(std::vector<cudaq::real>{probability}) {}
858+
REGISTER_KRAUS_CHANNEL()
859+
};
860+
600861
} // namespace cudaq

‎runtime/nvqir/cutensornet/simulator_cutensornet.cpp

+46
Original file line numberDiff line numberDiff line change
@@ -170,6 +170,52 @@ void SimulatorTensorNetBase::applyKrausChannel(
170170
}
171171
}
172172

173+
bool SimulatorTensorNetBase::isValidNoiseChannel(
174+
const cudaq::noise_model_type &type) const {
175+
switch (type) {
176+
case cudaq::noise_model_type::depolarization_channel:
177+
case cudaq::noise_model_type::bit_flip_channel:
178+
case cudaq::noise_model_type::phase_flip_channel:
179+
case cudaq::noise_model_type::x_error:
180+
case cudaq::noise_model_type::y_error:
181+
case cudaq::noise_model_type::z_error:
182+
case cudaq::noise_model_type::phase_damping:
183+
case cudaq::noise_model_type::pauli1:
184+
case cudaq::noise_model_type::pauli2:
185+
case cudaq::noise_model_type::depolarization1:
186+
case cudaq::noise_model_type::depolarization2:
187+
case cudaq::noise_model_type::unknown: // may be unitary, so return true
188+
return true;
189+
// These are explicitly non-unitary and unsupported
190+
case cudaq::noise_model_type::amplitude_damping_channel:
191+
case cudaq::noise_model_type::amplitude_damping:
192+
default:
193+
return false;
194+
}
195+
}
196+
197+
void SimulatorTensorNetBase::applyNoise(
198+
const cudaq::kraus_channel &channel,
199+
const std::vector<std::size_t> &targets) {
200+
LOG_API_TIME();
201+
// Do nothing if no execution context
202+
if (!executionContext)
203+
return;
204+
205+
// Do nothing if no noise model
206+
if (!executionContext->noiseModel)
207+
return;
208+
209+
// Apply all prior gates before applying noise.
210+
std::vector<int32_t> qubits{targets.begin(), targets.end()};
211+
cudaq::info(
212+
"[SimulatorTensorNetBase] Applying kraus channel {} on qubits: {}",
213+
cudaq::get_noise_model_type_name(channel.noise_type), qubits);
214+
215+
flushGateQueue();
216+
applyKrausChannel(qubits, channel);
217+
}
218+
173219
void SimulatorTensorNetBase::applyNoiseChannel(
174220
const std::string_view gateName, const std::vector<std::size_t> &controls,
175221
const std::vector<std::size_t> &targets,

‎runtime/nvqir/cutensornet/simulator_cutensornet.h

+6
Original file line numberDiff line numberDiff line change
@@ -36,6 +36,12 @@ class SimulatorTensorNetBase : public nvqir::CircuitSimulatorBase<double> {
3636
const std::vector<std::size_t> &targets,
3737
const std::vector<double> &params) override;
3838

39+
bool isValidNoiseChannel(const cudaq::noise_model_type &type) const override;
40+
41+
/// @brief Apply the given kraus_channel on the provided targets.
42+
void applyNoise(const cudaq::kraus_channel &channel,
43+
const std::vector<std::size_t> &targets) override;
44+
3945
// Override base calculateStateDim (we don't instantiate full state vector in
4046
// the tensornet backend). When the user want to retrieve the state vector, we
4147
// check if it is feasible to do so.

‎runtime/nvqir/stim/StimCircuitSimulator.cpp

+21-5
Original file line numberDiff line numberDiff line change
@@ -44,14 +44,30 @@ class StimCircuitSimulator : public nvqir::CircuitSimulatorBase<double> {
4444
isValidStimNoiseChannel(const kraus_channel &channel) const {
4545

4646
// Check the old way first
47-
if (channel.noise_type == cudaq::noise_model_type::bit_flip_channel)
47+
switch (channel.noise_type) {
48+
case cudaq::noise_model_type::bit_flip_channel:
49+
case cudaq::noise_model_type::x_error:
4850
return "X_ERROR";
49-
50-
if (channel.noise_type == cudaq::noise_model_type::phase_flip_channel)
51+
case cudaq::noise_model_type::y_error:
52+
return "Y_ERROR";
53+
case cudaq::noise_model_type::phase_flip_channel:
54+
case cudaq::noise_model_type::phase_damping:
55+
case cudaq::noise_model_type::z_error:
5156
return "Z_ERROR";
52-
53-
if (channel.noise_type == cudaq::noise_model_type::depolarization_channel)
57+
case cudaq::noise_model_type::depolarization_channel:
58+
case cudaq::noise_model_type::depolarization1:
5459
return "DEPOLARIZE1";
60+
case cudaq::noise_model_type::depolarization2:
61+
return "DEPOLARIZE2";
62+
case cudaq::noise_model_type::pauli1:
63+
return "PAULI_CHANNEL_1";
64+
case cudaq::noise_model_type::pauli2:
65+
return "PAULI_CHANNEL_2";
66+
case cudaq::noise_model_type::amplitude_damping_channel:
67+
case cudaq::noise_model_type::amplitude_damping:
68+
case cudaq::noise_model_type::unknown:
69+
return std::nullopt;
70+
}
5571

5672
return std::nullopt;
5773
}

0 commit comments

Comments
 (0)