This notebook serves as a compilation of the simulations which were run for the MMWU in Quantum Zero Sum Games paper.
If you wish to run the code yourself, please run pip install -r requirements.txt in order to install the required packages.
We also provide a html file which contains the code and outputs in lieu of running the code.
In [27]:
import pandas as pd
import numpy as np
import matplotlib
from matplotlib import pyplot as plt
import matplotlib.animation as animation
from mpl_toolkits import mplot3d
from scipy.integrate import odeint, solve_ivp
from scipy.stats import entropy
from scipy.linalg import sinm, cosm, logm, expm
import plotly.graph_objects as go
from plotly.subplots import make_subplots
from pathlib import Path
from sympy.physics.quantum import TensorProduct
from tqdm import tqdm
import cmath
import itertools
import warnings
from pylab import *
import matplotlib.animation as animation
from mpl_toolkits.mplot3d import Axes3D
from qutip import *
from algorithms import *
from helpers import *
%matplotlib inline
MMWU Simulations¶
In [2]:
def RunSimulation(payoff, basis1, basis2, basis3, basis4, s, discrete=False, discrete_num=5000, rep_num=500,
qre_num=500, alt=False, exponent=1/2, plot=True):
basis_data = GetBasisTransform(payoff=payoff, basis1=basis1, basis2=basis2, basis3=basis3, basis4=basis4, s=s)
print('\n')
print('Running discrete algorithm to obtain approximate Nash...')
print('\n')
discrete_data = RunParallelAlgoEpsDecrease(basis_data['R'],s=basis_data['s'], exponent=exponent, n=payoff.shape[0], m=payoff.shape[0], N=discrete_num, alt=alt)
discrete_traj = np.array([np.array([discrete_data['rho'][i], discrete_data['sig'][i]]) for i in range(len(discrete_data['rho']))])
print('Solving replicator system...')
print('\n')
rep_data = rep_trajectory(payoff=basis_data['R'], s = basis_data['s'], f=f_deriv, numstep=rep_num, n=payoff.shape[0], m=payoff.shape[0], plot=plot)
rep_traj = np.array([i.reshape(2,2,2) for i in rep_data['traj']])
print('Computing quantum relative entropies...')
print('\n')
if plot:
if discrete:
discrete_data_new = RunParallelAlgoEpsDecrease(basis_data['R'],s=basis_data['s'], exponent=exponent, n=payoff.shape[0], m=payoff.shape[0], N=qre_num)
PlotQRE(discrete_data_new, nash1=discrete_data['nash 1'], nash2=discrete_data['nash 2'], num=qre_num, discrete=discrete)
else:
PlotQRE(rep_data, nash1=discrete_data['nash 1'], nash2=discrete_data['nash 2'], num=qre_num, discrete=discrete)
print('Done!')
return discrete_traj, rep_traj
In [3]:
payoff = np.array([[1,-1],[-1,1]])
basis1 = np.matrix([[1, 0], [0, 1]])
basis2 = (1/np.sqrt(5))*np.matrix([[1, 2j], [2j, 1]])
basis3 = (1/np.sqrt(5))*np.matrix([[2, 1j], [1j, 2]])
In [4]:
sim_data_discrete1, sim_data_rep1 = RunSimulation(payoff, basis2, basis2, basis2, basis2, s=[0.35, 0.65, 0.65, 0.35], discrete=True, alt=False, discrete_num=50000, rep_num=500, exponent=1/4)
In [5]:
sig_data1 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete1[:,0]])
rho_data1 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete1[:,1]])
In [6]:
fig1 = go.Figure([go.Scatter(x = rho_data1[:,1], y=sig_data1[:,1],
mode='lines', line=dict(width=0.7, color='darkslateblue'),
name='Trajectories (Heads)'),
go.Scatter(x = [rho_data1[:,1][0]], y=[sig_data1[:,1][0]],
mode='markers',
name='Initial Condition', marker = dict(size=6, color='#440154', symbol = 'cross'), opacity=1),
])
# Edit layout
fig1.update_layout(title=r'$\mu = \log(1+ \frac{1}{t^{1/4}})$',
xaxis_title=r'$\sigma$',
yaxis_title=r'$\rho$',
xaxis_range=[0,1],
yaxis_range=[0,1],
legend_orientation='h',
legend=dict(y=-0.2),
font=dict(size=15),
width=500,
height=500,)
In [7]:
sim_data_discrete2, sim_data_rep2 = RunSimulation(payoff, basis2, basis2, basis2, basis2, s=[0.35, 0.65, 0.65, 0.35], discrete=True, alt=False, discrete_num=50000, rep_num=500, exponent=1/3)
In [8]:
sig_data2 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete2[:,0]])
rho_data2 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete2[:,1]])
In [9]:
fig2 = go.Figure([go.Scatter(x = rho_data2[:,1], y=sig_data2[:,1],
mode='lines', line=dict(width=0.7, color='darkslateblue'),
name='Trajectories (Heads)'),
go.Scatter(x = [rho_data2[:,1][0]], y=[sig_data2[:,1][0]],
mode='markers',
name='Initial Condition', marker = dict(size=6, color='#440154', symbol = 'cross'), opacity=1),
])
# Edit layout
fig2.update_layout(title=r'$\mu = \log(1+ \frac{1}{t^{1/4}})$',
xaxis_title=r'$\sigma$',
yaxis_title=r'$\rho$',
xaxis_range=[0,1],
yaxis_range=[0,1],
legend_orientation='h',
legend=dict(y=-0.2),
font=dict(size=15),
width=500,
height=500,)
In [10]:
sim_data_discrete3, sim_data_rep3 = RunSimulation(payoff, basis2, basis2, basis2, basis2, s=[0.35, 0.65, 0.65, 0.35], discrete=True, alt=False, discrete_num=50000, rep_num=500, exponent=1/2)
In [11]:
sig_data3 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete3[:,0]])
rho_data3 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete3[:,1]])
In [12]:
fig = go.Figure([go.Scatter(x = rho_data3[:,1], y=sig_data3[:,1],
mode='lines', line=dict(width=0.7, color='darkslateblue'),
name='Trajectories (Heads)'),
go.Scatter(x = [rho_data3[:,1][0]], y=[sig_data3[:,1][0]],
mode='markers',
name='Initial Condition', marker = dict(size=6, color='#440154', symbol = 'cross'), opacity=1),
])
# Edit layout
fig.update_layout(title=r'$\mu = \log(1+ \frac{1}{t^{1/2}})$',
xaxis_title=r'$\sigma$',
yaxis_title=r'$\rho$',
xaxis_range=[0,1],
yaxis_range=[0,1],
legend_orientation='h',
legend=dict(y=-0.2),
font=dict(size=15),
width=500,
height=500,)
In [13]:
sim_data_discrete4, sim_data_rep4 = RunSimulation(payoff, basis2, basis2, basis2, basis2, s=[0.35, 0.65, 0.65, 0.35], discrete=True, alt=False, discrete_num=50000, rep_num=500, exponent=2/3)
In [14]:
sig_data4 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete4[:,0]])
rho_data4 = np.array([np.linalg.eigvals(i).real for i in sim_data_discrete4[:,1]])
In [15]:
fig4 = go.Figure([go.Scatter(x = rho_data4[:,1], y=sig_data4[:,1],
mode='lines', line=dict(width=0.7, color='darkslateblue'),
name='Trajectories (Heads)'),
go.Scatter(x = [rho_data4[:,1][0]], y=[sig_data4[:,1][0]],
mode='markers',
name='Initial Condition', marker = dict(size=6, color='#440154', symbol = 'cross'), opacity=1),
])
# Edit layout
fig4.update_layout(title=r'$\mu = \log(1+ \frac{1}{t^{2/3}})$',
xaxis_title=r'$\sigma$',
yaxis_title=r'$\rho$',
xaxis_range=[0,1],
yaxis_range=[0,1],
legend_orientation='h',
legend=dict(y=-0.2),
font=dict(size=15),
width=500,
height=500,)
Constant of Motion¶
In [16]:
sim_data_discrete5, sim_data_rep5 = RunSimulation(payoff, basis3, basis2, basis3, basis2, s=[0.4, 0.6, 0.1, 0.9], discrete=False, alt=False, discrete_num=500000, rep_num=500)
Bloch Sphere¶
In [17]:
payoff_bloch1 = np.array([[1,-0.7],[-1,1.2]])
basis1 = np.matrix([[1, 0], [0, 1]])
basis2 = (1/np.sqrt(2))*np.matrix([[1, 1j], [1j, 1]])
basis3 = (1/np.sqrt(5))*np.matrix([[1, 2j], [2j, 1]])
In [18]:
sim_data_discrete_bloch1, sim_data_rep_bloch1 = RunSimulation(payoff_bloch1, basis3, basis2, basis3, basis2, s=[0.2, 0.8, 0.7, 0.3], discrete_num=50000, discrete=False, alt=False, plot=False)
In [19]:
states1 = [Qobj(sim_data_rep_bloch1[i][0]) for i in range(500)]
states2 = [Qobj(sim_data_rep_bloch1[i][1]) for i in range(500)]
# AnimateBloch(states1, states2, duration=0.05, save_all=False, filename='bloch_rep_1', length=500)
In [20]:
payoff_bloch2 = np.array([[1,0], [-1, 1]])
basis1 = (1/np.sqrt(2))*np.matrix([[1, 1], [-1, 1]])
basis2 = (1/np.sqrt(5))*np.matrix([[1, 2], [-2, 1]])
In [21]:
sim_data_discrete_bloch2, sim_data_rep_bloch2 = RunSimulation(payoff_bloch2, basis1, basis2, basis1, basis2, s=[0.45, 0.55, 0.4, 0.6], discrete_num=50000, discrete=False, alt=False, rep_num=1000, qre_num=500, plot=False)
In [22]:
states3 = [Qobj(sim_data_rep_bloch2[i][0]) for i in range(1000)]
states4 = [Qobj(sim_data_rep_bloch2[i][1]) for i in range(1000)]
# AnimateBloch(states3, states4, duration=0.01, save_all=True, filename='bloch_rep_2', length=500)
In [23]:
PrintBloch(states3, states4, timestep=0, save=False)
In [24]:
PrintBloch(states3, states4, timestep=49, save=False)
In [25]:
PrintBloch(states3, states4, timestep=369, save=False)
In [26]:
PrintBloch(states3, states4, timestep=999, save=False)
Bloch Sphere animation showing Poincare recurrence in the quantum setting.
