Skip to content
Permalink
7bd854fbf6
Switch branches/tags

Name already in use

A tag already exists with the provided branch name. Many Git commands accept both tag and branch names, so creating this branch may cause unexpected behavior. Are you sure you want to create this branch?

B92Finite/utils.cpp

Go to file
 
 
Cannot retrieve contributors at this time
391 lines (277 sloc) 8.88 KB
Raw Blame
#include "utils.hpp"
namespace Algebra{
using namespace std;
void eigenValuestoVectorRef(QMat m, vector<complex<double>>& v)
{
QMat e = m.eigenValues();
for(int i=0; i<e.rows(); ++i){
v.push_back(e.getElement(i,0));
}
}
map<string, double> random_angles()
{
random_device rd;
mt19937 gen(rd());
uniform_real_distribution<> twopi(0.0, 2 * PI);
uniform_real_distribution<> halfpi(0.0, PI/2);
map<string, double> angles;
angles["phi"] = halfpi(gen);
angles["psi"] = twopi(gen);
angles["chi"] = twopi(gen);
return angles;
}
map<string, double> near_ident_random_angles()
{
random_device rd;
mt19937 gen(rd());
double t = T;
uniform_real_distribution<> dist(-t, t);
map<string, double> angles;
angles["phi"] = dist(gen);
angles["psi"] = dist(gen);
angles["chi"] = dist(gen);
return angles;
}
QMat random_unitary_matrix(int d) // dimension
{
vector<QMat> E;
for(int n=0; n<d-1; n++){
vector<QMat> compRots;
for(int i=n; i>=0; i--){
QMat cr(d,true);
map<string,double> angles = near_ident_random_angles();
complex<double> ii(0, angles["psi"]);
cr.setElement(i,i,cos(angles["phi"]) * exp(ii));
complex<double> jj(0,-1*angles["psi"]);
cr.setElement(n+1,n+1,cos(angles["phi"]) * exp(jj));
if(i == 0){
complex<double> ij(0,angles["chi"]);
cr.setElement(i,n+1,sin(angles["phi"]) * exp(ij));
complex<double> ji(0,-1*angles["chi"]);
cr.setElement(n+1,i, -1 * sin(angles["phi"]) * exp(ji));
}else{
complex<double> ij(0,0);
cr.setElement(i,n+1, sin(angles["phi"]) * exp(ij));
complex<double> ji(0,0);
cr.setElement(n+1,i, -1 * sin(angles["phi"]) * exp(ji));
}
compRots.push_back(cr);
}
QMat En = compRots[0];
for(int m=1; m<compRots.size(); m++){
En*=compRots[m];
}
E.push_back(En);
}
QMat U = E[0];
for(int m=1; m<E.size(); m++){
U*=E[m];
}
return U;
}
map<string, double> noiseStatsAtoB(QMat& Uf, vector<Qubit>& iSet, vector<Qubit>& jSet)
{
map<string, double> noises;
QMat zeroM(16,1);
complex<double> one0(1.0,0.0);
zeroM.setElement(0,0,one0);
for(auto& i : iSet){
QMat i0 = i.second.tensor(zeroM);
QMat i0star = i0.conjugateTranspose();
QMat Ufstar = Uf.conjugateTranspose();
QMat right = Uf * i0 * i0star * Ufstar;
for(auto& j : jSet){
QMat jstar = j.second.conjugateTranspose();
QMat jjstar = j.second*jstar;
QMat Ie(16, true);
QMat left = jjstar.tensor(Ie);
QMat P = left * right;
string qtq = "P" + i.first + j.first;
noises[qtq] = P.trace().real();
}
}
return noises;
}
map<string, double> noiseStatsAtoBreftoA(QMat& Uf, QMat& Ur, vector<Qubit>& iSet, vector<Qubit>& jSet)
{
map<string, double> noises;
QMat zeroM(16,1);
complex<double> one0(1.0,0.0);
zeroM.setElement(0,0,one0);
for(auto& i : iSet){
QMat i0 = i.second.tensor(zeroM);
QMat i0star = i0.conjugateTranspose();
QMat Ufstar = Uf.conjugateTranspose();
QMat Urstar = Ur.conjugateTranspose();
QMat right = Ur* Uf * i0 * i0star * Ufstar * Urstar;
for(auto& j : jSet){
QMat jstar = j.second.conjugateTranspose();
QMat jjstar = j.second*jstar;
QMat Ie(16, true);
QMat left = jjstar.tensor(Ie);
QMat P = left * right;
string qtq = "P" + i.first + "R" + j.first;
noises[qtq] = P.trace().real();
}
}
return noises;
}
map<string, double> noiseStatsAtoBprojtoA(QMat& Uf, QMat& Ur, vector<Qubit>& iSet, vector<Qubit>& jSet, vector<Qubit>& kSet)
{
map<string, double> noises;
QMat zeroM(16,1);
complex<double> one0(1.0,0.0);
zeroM.setElement(0,0,one0);
for(auto& i : iSet){
QMat i0 = i.second.tensor(zeroM);
QMat i0star = i0.conjugateTranspose();
QMat Ufstar = Uf.conjugateTranspose();
QMat Urstar = Ur.conjugateTranspose();
QMat right = Uf * i0 * i0star * Ufstar;
for(auto& j : jSet){
QMat jstar = j.second.conjugateTranspose();
QMat jjstar = j.second*jstar;
QMat Ie(16, true);
QMat left = jjstar.tensor(Ie);
QMat P = left * right;
complex<double> trs(1.0 / P.trace().real(), 0.0);
P*=left;
P*=trs;
for(auto& k: kSet){
QMat kstar = k.second.conjugateTranspose();
QMat kkstar = k.second*kstar;
QMat kleft = kkstar.tensor(Ie);
QMat F = kleft * Ur * P * Urstar;
string qtqtq = "P" + i.first + j.first + k.first;
noises[qtqtq] = F.trace().real();
}
}
}
return noises;
}
map<string, double> true_value(QMat& Uf, QMat& Ur, Qubit& zero, Qubit& one)
{
complex<double> one0(1.0,0.0);
QMat first(32,1);
first.setElement(0,0,one0);
QMat second(32,1);
second.setElement(16,0,one0);
QMat Uffirst = Uf*first;
QMat Ufsecond = Uf*second;
QMat e0(16,1);
QMat e1(16,1);
QMat e2(16,1);
QMat e3(16,1);
for(int i=0; i<16; i++){
e0.setElement(i,0,Uffirst.getElement(i,0));
e2.setElement(i,0,Ufsecond.getElement(i,0));
}
for(int i=16; i<32; i++){
e1.setElement(i-16,0,Uffirst.getElement(i,0));
e3.setElement(i-16,0,Ufsecond.getElement(i,0));
}
QMat e00 = Ur * zero.second.tensor(e0);
QMat e11 = Ur * one.second.tensor(e1);
QMat e02 = Ur * zero.second.tensor(e2);
QMat e13 = Ur * one.second.tensor(e3);
QMat e000(16,1);
QMat e110(16,1);
QMat e020(16,1);
QMat e130(16,1);
QMat e001(16,1);
QMat e111(16,1);
QMat e021(16,1);
QMat e131(16,1);
for(int i=0; i<16; i++){
e000.setElement(i,0,e00.getElement(i,0));
e110.setElement(i,0,e11.getElement(i,0));
e020.setElement(i,0,e02.getElement(i,0));
e130.setElement(i,0,e13.getElement(i,0));
}
for(int i=16; i<32; i++){
e001.setElement(i-16,0,e00.getElement(i,0));
e111.setElement(i-16,0,e11.getElement(i,0));
e021.setElement(i-16,0,e02.getElement(i,0));
e131.setElement(i-16,0,e13.getElement(i,0));
}
map <string,double> values;
complex<double> E1 = e000.inner_product(e131);
values["E1"] = E1.real();
complex<double> E2 = e111.inner_product(e020);
values["E2"] = E2.real();
complex<double> E3 = e001.inner_product(e130);
values["E3"] = E3.real();
complex<double> E4 = e110.inner_product(e021);
values["E4"] = E4. real();
return values;
}
};
int maain()
{
using namespace Algebra;
using namespace std;
complex<double> one0(1.0,0.0);
complex<double> oneoRoot2(1/sqrt(2),0.0);
complex<double> noneoRoot2(-1.0/sqrt(2),0.0);
complex<double> ioRoot2(0.0,1.0/sqrt(2));
complex<double> nioRoot2(0.0,-1.0/sqrt(2));
QMat zeroM(2,1);
zeroM.setElement(0,0,one0);
Qubit zero("0", zeroM);
QMat oneM(2,1);
oneM.setElement(1,0,one0);
Qubit one("1", oneM);
QMat plusM(2,1);
plusM.setElement(0,0,oneoRoot2);
plusM.setElement(1,0,oneoRoot2);
Qubit plus("P", plusM);
QMat minusM(2,1);
minusM.setElement(0,0,oneoRoot2);
minusM.setElement(1,0,noneoRoot2);
Qubit minus("M", minusM);
QMat zeroyM(2,1);
zeroyM.setElement(0,0,oneoRoot2);
zeroyM.setElement(1,0,ioRoot2);
Qubit zeroy("ZY", zeroyM);
QMat oneyM(2,1);
oneyM.setElement(0,0,oneoRoot2);
oneyM.setElement(1,0,nioRoot2);
Qubit oney("OY", oneyM);
vector<Qubit> iSet = {zero, one, plus};
vector<Qubit> jSet = {zero, one};
vector<Qubit> kSet = {zero, one, plus, minus};
QMat Uf = random_unitary_matrix(32);
QMat Ur = random_unitary_matrix(32);
map<string, double> AtoBStats = noiseStatsAtoB(Uf, iSet, jSet);
map<string, double> AtoBreftoAStats = noiseStatsAtoBreftoA(Uf, Ur, iSet, kSet);
map<string, double> AtoBprojtoAStats = noiseStatsAtoBprojtoA(Uf,Ur, iSet, jSet, kSet);
map<string, double> values = true_value(Uf, Ur, zero, one);
cout << "Uf:" << endl;
Uf.serialize();
cout << endl;
cout << "Ur:" << endl;
Ur.serialize();
cout << endl;
cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
cout << "A -> B stats:" << endl;
for (auto& ns : AtoBStats){
cout << ns.first << " = " << ns.second << endl;
}
cout << endl;
cout << "A -> B reflects -> A stats:" << endl;
for (auto& ns : AtoBreftoAStats){
cout << ns.first << " = " << ns.second << endl;
}
cout << endl;
cout << "A -> B projection -> A stats:" << endl;
for (auto& ns : AtoBprojtoAStats){
cout << ns.first << " = " << ns.second << endl;
}
cout << endl;
cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
cout << "True value:" << endl;
for (auto& ns : values){
cout << ns.first << " = " << ns.second << endl;
}
return 0;
}