// quantum computing emulation: Super Dense Coding with phase shift error
//
// Usage: SDC [-rdbE] [rpA [rpB]]
//
//  rp = phase shift error relative values for Alice and Bob, default = 1 (no error)
//  -r  use random global phase factor
//  -d  extra debug output, -dd for even more output
//  -b  print state indices in binary
//  -E  show entanglement
//
// R. Perry, June 2020
//
#include <iostream>
#include <cstdio>  // sprintf()
#include <cstdlib> // atof()
#include "qce.h"

using namespace qce; using std::cout; // using namespace std;

#define Real double
#define Complex std::complex<Real>
#define State state<Complex,Real>

// debug: display transformation result
//
void disp( const State &q, const char *what, int c, int k, print_option b, int E)
{
  char msg[20];

  if( c >= 0) sprintf( msg, "%s%i%i", what, c, k);
  else if( k >= 0) sprintf( msg, "%s%i", what, k);
  else sprintf( msg, "%s", what);

  q.print( msg, b); if( E) q.print( "E", PRINT_TANGLE);
}

int main( int argc, char *argv[])
{
  int rphase = 0, debug = 0, E = 0; unsigned int seed = SRAND();

  print_option b = PRINT_NONE;

 // command-line options
 //
  if( argc > 1 && argv[1][0] == '-') 
  {
    for( int i = 1; argv[1][i]; ++i) switch( argv[1][i])
    {
      case 'r': ++rphase; break;
      case 'd': ++debug; break;
      case 'b': b = PRINT_BITS; break;
      case 'E': ++E; break;
      default: error( "SDC: bad option");
    }
    --argc, ++argv;
  }

 // phase shift relative value
 //
  Real rpA = 1, rpB = 1;

  if( argc > 1) { rpA = atof(argv[1]); if( argc > 2) rpB = atof(argv[2]); }

  if( rpA < 0 || rpA > 2 || rpB < 0 || rpB > 2) error( "SDC: bad rp arg");

  State q(2);

  cout << "SDC: seed = " << seed << ", rpA = " << rpA << ", rpB = " << rpB << "\n";

 // encode value 0 ... 3
 //
  for( unsigned int i = 0; i < 4; ++i)
  {
    cout << i << ":\n"; 

   //----------------------------------------------------------------------
   // Alice+Bob: initialize, H, CX
   //
    q.init( 0, rphase); if( debug > 1) disp( q, "q", -1, -1, b, E);
   //
   // apply Hrp to qubit 1
   //
    q.Hrp( 1, rpA); if( debug > 1) disp( q, "Hrp", -1, 1, b, E);
   //
   // apply CX to qubit 0 controlled by qubit 1
   //
    q.CX( 1, 0); if( debug > 1) disp( q, "CX", 1, 0, b, E);

   if( debug == 1) disp( q, "A", -1, -1, b, E);

   //----------------------------------------------------------------------
   // Alice: encode i using X and Z transformations on qubit 1
   //
    if( i & 1)
    {
      q.X( 1); if( debug) disp( q, "X", -1, 1, b, E);
    }
   //
    if( i & 2)
    {
      q.Z( 1); if( debug) disp( q, "Z", -1, 1, b, E);
    }
   //
   // send qubit 1 to Bob
    
   //----------------------------------------------------------------------
   // Bob: CX, H, measure
   //
   // apply CX to qubit 0 controlled by qubit 1
   //
    cout << "Bob:\n";
   //
    q.CX( 1, 0); if( debug > 1) disp( q, "CX", 1, 0, b, E);

   // apply H to qubit 1
   //
    q.Hrp( 1, rpB); if( debug > 1) disp( q, "Hrp", -1, 1, b, E);
 
    if( debug == 1) disp( q, "B", -1, -1, b, E);

   // measure: check max probability
   //
    unsigned int s = 0; Real p = abs2(q.a[0]), t;

    for( unsigned int j = 1; j < 4; ++j) { t = abs2(q.a[j]); if( t > p) { p = t; s = j; } }

    cout << " i = " << i << " , P(" << s << ") = " << p << "\n";
   //
   //----------------------------------------------------------------------
   //
  }

  return 0;
}