Notebook for creating the figures in "Quantum Computing and Tensor Networks for Laminate Design: A Novel Approach to Stacking Sequence Retrieval"¶
Disclaimer: This repository contains experimental software and is published for the sole purpose of giving additional background details on the respective publication.
In [1]:
import numpy as np
from numpy import pi
import h5py
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib.colors import hsv_to_rgb
from os import listdir
from os.path import isfile,join
from copy import deepcopy
from scipy.optimize import curve_fit
In [2]:
# matplotlib settings
mpl.rcParams['font.sans-serif'] = 'Arial'
mpl.rcParams['font.family'] = 'Arial'
mpl.rcParams['axes.linewidth'] = 0.6
mpl.rcParams['xtick.major.size'] = 2.5
mpl.rcParams['xtick.major.width'] = 0.6
mpl.rcParams['xtick.minor.size'] = 1.5
mpl.rcParams['xtick.minor.width'] = 0.6
mpl.rcParams['ytick.major.size'] = 2.5
mpl.rcParams['ytick.major.width'] = 0.6
mpl.rcParams['ytick.minor.size'] = 1.5
mpl.rcParams['ytick.minor.width'] = 0.6
mpl.rcParams.update({'font.size': 10})
plt.rcParams['svg.fonttype'] = 'none'
Function definitions¶
In [3]:
# function for calculation of lamiantion parameters
def lamination_parameters(angles, symmetric=True, num_plies_odd=False, parameters='AD', deg=False):
# convert angles, if in degrees
if deg:
angles = pi / 180 * np.array(angles)
# mirror angles if ply is symmetric
if symmetric:
if num_plies_odd:
angles = np.concatenate((angles[::-1], angles[1:]))
else:
angles = np.concatenate((angles[::-1], angles))
num_plies = len(angles)
dz = 1 / num_plies
boundaries = -1 / 2 + np.arange(num_plies + 1) * dz
lps = np.array([])
for par in parameters:
if par.upper() == 'A':
weights = np.full((num_plies,), dz)
fac = 1
elif par.upper() == 'B':
weights = boundaries[1:] ** 2 - boundaries[:-1] ** 2
fac = 2
elif par.upper() == 'D':
weights = boundaries[1:] ** 3 - boundaries[:-1] ** 3
fac = 4
else:
raise ValueError(f'Argument "parameters" should only contain A,B and D, not {par}')
lp = fac * np.array([np.sum(weights * np.cos(2 * angles)),
np.sum(weights * np.sin(2 * angles)),
np.sum(weights * np.cos(4 * angles)),
np.sum(weights * np.sin(4 * angles))])
lps = np.concatenate((lps, lp))
return lps
# function that visualizes a stacking sequence
def plot_laminate(stack,angles,zero_indexed=True,save_fig=False, plot_target=False, plot_disor_constr_violations=False):
num_plies = len(stack)
num_angles = len(angles)
color_list = [hsv_to_rgb((n / num_angles, 1, 1)) for n in range(num_angles)]
fig, ax = plt.subplots(figsize=(8, 1.75))
ax.set_frame_on(False)
ax.get_xaxis().tick_bottom()
ax.axes.get_yaxis().set_visible(False)
ax.set_xlim(0, num_plies)
if not isinstance(stack,np.ndarray):
stack = np.array(stack)
if not zero_indexed:
stack = stack-1
for k in range(num_angles):
x = np.where(stack == k)[0]
ax.bar(x, 1., width=1., color=color_list[k], label=f'{angles[k]}°', align='edge')
if num_angles <= 4:
ncol = 1
else:
ncol = 2
plt.legend(bbox_to_anchor=(1.01, 1), loc='upper left', borderaxespad=0, ncol=ncol)
if plot_target:
num_t_angles = len(target_angles_indices)
ymax = 1.1
for k in range(num_angles):
x = np.where(target_angles_indices == k)[0]
ax.bar(x * num_plies / num_t_angles, ymax, width=num_plies / num_t_angles, bottom=1.025,
color=color_list[k],
label=f'{angles[k]}°', align='edge')
else:
ymax = 1.
if plot_disor_constr_violations:
ymin = -0.1
ax.bar([0.], ymin, width=num_plies, bottom=-0.025, color=(0.8, 0.8, 0.8), align='edge')
ax.bar(disor_constr_violations, ymin, width=1., bottom=-0.025, color='red', align='edge')
else:
ymin = 0.
ax.set_ylim(ymin, ymax)
return fig, ax
# wrapper for lamination parameters for symmetric stack with angles in degrees
def lp_short(angles):
return lamination_parameters(angles,symmetric=True,num_plies_odd=False,parameters='AD',deg=True)
# calculate LP a stack in terms of angle indices, from angles
def stack_to_lp(stack,angles):
return lamination_parameters(angles[stack-1],symmetric=True,num_plies_odd=False,parameters='AD',deg=True)
# calculate loss function between lamination parameters
def lp_to_loss(lp,target):
return np.sum((lp-target)**2)
In [4]:
# Functions for accessing data from the HDF5 files
# generate folder path for specific configurations
def folderpath(filelocation,constraint,num_sweeps,max_bond_dim,sweep_direction):
filelocation = filelocation.replace('\\','/')
if filelocation[-1] != '/':
filelocation += '/'
if constraint:
filelocation += 'with_constraint/'
else:
filelocation += 'without_constraint/'
filelocation += f'num_sweeps_{num_sweeps}/max_bond_dim_{max_bond_dim}/sweep_{sweep_direction}/'
return filelocation
# get all the filenames in the folder
def get_files_in_folder(filelocation,filename,
constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None,return_full_path=False):
if (constraint is None
or num_sweeps is None
or max_bond_dim is None
or sweep_direction is None):
filelocation = filelocation.replace('\\','/')
if filelocation[-1] != '/':
filelocation += '/'
else:
filelocation = folderpath(filelocation,constraint,num_sweeps,max_bond_dim,sweep_direction)
files = [f for f in listdir(filelocation)
if isfile(join(filelocation,f)) and f.endswith('.hdf5') and f.startswith(filename+'_sample')]
if return_full_path:
files = [join(filelocation,f) for f in files]
return files
# get a file from the folder using the sample index
def get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None,return_full_path=False):
if (constraint is None
or num_sweeps is None
or max_bond_dim is None
or sweep_direction is None):
filelocation = filelocation.replace('\\','/')
if filelocation[-1] != '/':
filelocation += '/'
else:
filelocation = folderpath(filelocation,constraint,num_sweeps,max_bond_dim,sweep_direction)
filepath = filelocation + filename + f'_sample_{sample_idx:04}.hdf5'
if return_full_path:
return join(filelocation,filepath)
return filepath
# create a dictionary with the most important data from the HDF5 file
def file_read_main_results(filelocation,filename,sample_idx,
constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None):
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
res_dict = {}
with h5py.File(filepath,'r') as file:
res_dict['target_parameters'] = np.array(file['properties/target_parameters'])
res_dict['target_stack'] = np.array(file['properties/target_stack'])
res_dict['energies'] = np.array(file['data/energies'])
res_dict['elapsed_time_sweeps'] = np.array(file['data/elapsed_time_sweeps'])
res_dict['time_stamps'] = np.array(file['data/time_stamps'])
res_dict['results_loss'] = np.array(file['results/loss'])
res_dict['results_rmse'] = np.array(file['results/rmse'])
res_dict['results_stack'] = np.array(file['results/stack'])
res_dict['results_lamination_parameters'] = np.array(file['results/lamination_parameters'])
res_dict['results_constriant_violations'] = np.array(file['results/constraint_violations'])
res_dict['psi_samples'] = {}
for n in file['data/psi_samples'].keys():
res_dict['psi_samples'][n] = {}
for nds,ds in file['data/psi_samples'][n].items():
res_dict['psi_samples'][n][nds] = np.array(ds,dtype=int)
return res_dict
# get the loss function evaluations from the results of the optimizations
def get_losses_for_samples(filelocation,filename,sample_indices,constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None):
filepath = get_file_from_sample_idx(filelocation,filename,sample_indices[0],
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
num_tries = file['properties'].attrs['num_tries']
losses = np.empty((len(sample_indices),num_tries),dtype = float)
for i,sample_idx in enumerate(sample_indices):
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
losses[i,:] = np.array(file['results/loss'])
return losses
# get the root mean square errors from the results of the optimizations
def get_rmse_for_samples(filelocation,filename,sample_indices,constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None):
filepath = get_file_from_sample_idx(filelocation,filename,sample_indices[0],
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
num_tries = file['properties'].attrs['num_tries']
rmses = np.empty((len(sample_indices),num_tries),dtype = float)
for i,sample_idx in enumerate(sample_indices):
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
rmses[i,:] = np.array(file['results/rmse'])
return rmses
File location and name¶
In [5]:
filelocation = 'C:/Users/arnewulff/Documents/ExperimentalData/stacking_retrieval/constant_bd/'
filename = 'constant_bondsize'
Check constraint violations¶
In [6]:
# Function for counting constraint violations:
constr_matrix = np.array([[0,0,1,0],
[0,0,0,1],
[1,0,0,0],
[0,1,0,0]])
def count_constraint_violations(stack,zero_indexed=False):
if not zero_indexed:
stack = stack-1
viols = 0
for p1,p2 in zip(stack[:-1],stack[1:]):
viols += constr_matrix[p1,p2]
return viols
num_sweeps = 69
sample_indices = np.arange(50)
bond_dims = [2,4,8,16,32]
sweep_sequences = ['L','LR','R']
# counting constraint violations
# without enforcement of the constraint
print('**Without** enforcement of the constraint:')
constraint = False
num_constr_viols = 0
for r,sd in enumerate(sweep_sequences):
for c,bd in enumerate(bond_dims):
for sample_idx in sample_indices:
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,constraint=constraint,
num_sweeps=num_sweeps,max_bond_dim=bd,sweep_direction=sd)
with h5py.File(filepath,'r') as file:
stacks = np.array(file['results/stack'],dtype=int)
constraint_violations = np.apply_along_axis(count_constraint_violations,0,stacks,zero_indexed=False)
if np.sum(constraint_violations) > 0:
num_constr_viols += 1
print(f'Number of samples with constraint violations: {num_constr_viols}\n')
# with enforcement of the constraint
print('**With** enforcement of the constraint:')
constraint = True
num_constr_viols = 0
for r,sd in enumerate(sweep_sequences):
for c,bd in enumerate(bond_dims):
for sample_idx in sample_indices:
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,constraint=constraint,
num_sweeps=num_sweeps,max_bond_dim=bd,sweep_direction=sd)
with h5py.File(filepath,'r') as file:
stacks = np.array(file['results/stack'],dtype=int)
constraint_violations = np.apply_along_axis(count_constraint_violations,0,stacks,zero_indexed=False)
if np.sum(constraint_violations) > 0:
num_constr_viols += 1
print(f'Number of samples with constraint violations: {num_constr_viols}\n')
**Without** enforcement of the constraint: Number of samples with constraint violations: 733 **With** enforcement of the constraint: Number of samples with constraint violations: 0
Figure of all measured samples¶
In [7]:
# Function for plotting the best, worst and and average for each configuration and
def plot_measured_points(ax,measured_points,color='red',alpha=1.,lw=1,ms=1):
"""Plot points with average, best and worst.
ax: axis to plot on
measured_points: array with shape (num_samples,num_tries)
"""
best_points = measured_points.min(axis=1)
worst_points = measured_points.max(axis=1)
averages = measured_points.mean(axis=1)
x = np.arange(measured_points.shape[0])
for i,(bp,wp) in enumerate(zip(best_points,worst_points)):
ax.plot([i,i],[bp,wp],'-',color=color,alpha=alpha,lw=lw)
ax.plot(x,worst_points,'.',color=color,alpha=alpha,mfc='white',ms=ms,mew=lw)
ax.plot(x,best_points,'.',color=color,alpha=alpha,mfc='white',ms=ms,mew=lw)
ax.plot(x,averages,'.',color=color,alpha=alpha,ms=ms,mew=lw)
In [8]:
# Plot figure of all samples
cm=1/2.54
# compare bonddims
num_sweeps = 69
sample_indices = np.arange(50)
bond_dims = [2,4,8,16,32]
sweep_sequences = ['L','LR','R']
# get sorted sample indices by sorting for 'L', highest bd without constr.
measured_points = get_rmse_for_samples(
filelocation,filename,sample_indices,constraint=False,
num_sweeps=num_sweeps,max_bond_dim=bond_dims[-1],sweep_direction=sweep_sequences[0]
)
sorted_indices = np.argsort(measured_points.mean(axis=1))
sorted_sample_indices = sample_indices[sorted_indices]
# create plot
fig,axs = plt.subplots(len(sweep_sequences),len(bond_dims),sharex=True,sharey='row',figsize=(27*cm,18*cm))
fig.subplots_adjust(hspace=0.1,wspace=0)
for r,sd in enumerate(sweep_sequences):
axs[r,0].set_ylabel('RMSE')
for c,bd in enumerate(bond_dims):
measured_points = get_rmse_for_samples(
filelocation,filename,sorted_sample_indices,constraint=False,
num_sweeps=num_sweeps,max_bond_dim=bd,sweep_direction=sd
)
measured_points_constr = get_rmse_for_samples(
filelocation,filename,sorted_sample_indices,constraint=True,
num_sweeps=num_sweeps,max_bond_dim=bd,sweep_direction=sd
)
plot_measured_points(axs[r,c],measured_points_constr,color=(0.3,0.3,1.),alpha=1.,ms=3.2,lw=0.6)
plot_measured_points(axs[r,c],measured_points,color='red',alpha=1.,ms=3.2,lw=0.6)
axs[r,c].grid(axis='y',lw=0.6)
# top titles
l,b,w,h = axs[0,0].get_position().bounds
# tit.set_position((tit.get_position()[0],b+h+0.05))
fig.text(l,b+h+0.005,' Bond-dim.:',va='bottom',ha='left',weight='bold')
for c,bd in enumerate(bond_dims):
l,b,w,h = axs[0,c].get_position().bounds
if c == 0:
fig.text(l+w/2+0.005,b+h+0.005,f'{bd}',va='bottom',ha='center')
continue
fig.text(l+w/2,b+h+0.005,f'{bd}',va='bottom',ha='center')
# left titles
l,b,w,h = axs[0,-1].get_position().bounds
fig.text(l+w-0.038,b+h-0.015,f'Sweep-dir.:',va='top',ha='right',weight='bold')
fig.text(l+w-0.008,b+h-0.015,f'Left',va='top',ha='right')
l,b,w,h = axs[1,-1].get_position().bounds
fig.text(l+w-0.008,b+h-0.015,f'Left-Right',va='top',ha='right')
l,b,w,h = axs[2,-1].get_position().bounds
fig.text(l+w-0.008,b+h-0.015,f'Right',va='top',ha='right')
axs[2,0].set_yticks([0.,0.05,0.1,0.15])
axs[1,0].set_ylim((axs[1,0].get_ylim()[0],0.15))
axs[0,0].set_ylim((axs[0,0].get_ylim()[0],0.06))
axs[-1,0].set_xticks([25-1,50-1],[25,50])
axs[-1,0].set_xlim((-1.5,50.5))
# create legend
l,b,w,h = axs[0,-1].get_position().bounds
ax_legend = fig.add_axes([l+w,b,w/2,h])
ax_legend.set_xlim(0.,1.)
ax_legend.set_ylim(0.,1.05)
ax_legend.plot([0.1,0.1],np.array([0.5,0.9])-0.02,'-',color='grey',lw=0.6)
ax_legend.plot([0.1],[0.7-0.02],'.',color='grey',mew=0.6)
ax_legend.plot([0.1,0.1],np.array([0.5,0.9])-0.02,'.',color='grey',mfc='white',mew=0.6)
ax_legend.text(0.2,0.48,'best',va='center',ha='left')
ax_legend.text(0.2,0.68,'mean',va='center',ha='left')
ax_legend.text(0.2,0.88,'worst',va='center',ha='left')
ax_legend.plot([0.1,0.1],np.array([0.3,0.4])-0.02,'-',color='red',lw=0.6)
ax_legend.plot([0.1],[0.35-0.02],'.',color='red',mew=0.6)
ax_legend.plot([0.1,0.1],np.array([0.3,0.4])-0.02,'.',color='red',mfc='white',mew=0.6)
ax_legend.text(0.2,0.35-0.02,'w/o constraint',va='center',ha='left')
ax_legend.plot([0.1,0.1],np.array([0.1,0.2])-0.02,'-',color=(0.3,0.3,1.),lw=0.6)
ax_legend.plot([0.1],[0.15-0.02],'.',color=(0.3,0.3,1.),mew=0.6)
ax_legend.plot([0.1,0.1],np.array([0.1,0.2])-0.02,'.',color=(0.3,0.3,1.),mfc='white',mew=0.6)
ax_legend.text(0.2,0.15-0.02,'w/ constraint',va='center',ha='left')
ax_legend.text(0.05,1.00,'Legend:',va='center',ha='left',weight='bold')
ax_legend.axis('off')
# plt.savefig('...',bbox_inches='tight',dpi=500)
Out[8]:
(0.0, 1.0, 0.0, 1.05)
Figure for the average energy expectation values over the course of the optimization¶
In [9]:
# Function for getting the average energy expectation values
def get_mean_energy_traces(filelocation,filename,sample_indices,constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None,add_initial_energies=False):
filepath = get_file_from_sample_idx(filelocation,filename,sample_indices[0],
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
num_tries = file['properties'].attrs['num_tries']
if num_sweeps is None:
num_sweeps = file['properties'].attrs['num_sweeps']
if add_initial_energies:
num_sweeps_tot = num_sweeps + 1
fi = 1
else:
num_sweeps_tot = num_sweeps
fi = 0
energies = np.empty((len(sample_indices),num_sweeps_tot,num_tries))
for i,sample_idx in enumerate(sample_indices):
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
energies[i,fi:,:] = np.array(file['data/energies'])
if add_initial_energies:
for ti,t in enumerate(file['psi0']):
energies[i,0,ti] = file[f'psi0/{t}'].attrs['energy']
return energies.transpose((1,0,2)).reshape(num_sweeps_tot,-1).mean(axis=1)
In [10]:
cm = 1/2.54
num_sweeps=69
sweep_indices = np.arange(50)
fig,axs = plt.subplots(1,3,figsize=(16*cm,8*cm),sharey=True)
fig.subplots_adjust(wspace=0)
for c,sd in enumerate(['L','LR','R']):
axs[c].set_yscale('log')
axs[c].grid(lw=0.6)
axs[c].set_xlim(-2,82)
axs[c].set_xticks([0,20,40,60])
for bdi,bd in enumerate(bond_dims):
ens = get_mean_energy_traces(filelocation,filename,sweep_indices,constraint=True,num_sweeps=num_sweeps,max_bond_dim=bd,
sweep_direction=sd,add_initial_energies=True)
sidx = 69-int(np.log2(bd))
if c == 2 and bdi==0:
lab = 'w/ constraint'
else:
lab = None
axs[c].plot(np.arange(sidx),ens[:sidx],'-',color=(0.3,0.3,1.),label=lab,lw=0.6)
axs[c].plot(np.arange(sidx-1,70),ens[sidx-1:],':',color=(0.3,0.3,1.),lw=0.6)
axs[c].text(74,ens[-1],bd,va='center',ha='left',color=(0.3,0.3,1.),size=8)
for bdi,bd in enumerate(bond_dims):
ens = get_mean_energy_traces(filelocation,filename,sweep_indices,constraint=False,num_sweeps=num_sweeps,max_bond_dim=bd,
sweep_direction=sd,add_initial_energies=True)
sidx = 69-int(np.log2(bd))
if c == 2 and bdi==0:
lab = 'w/o constraint'
else:
lab = None
axs[c].plot(np.arange(sidx),ens[:sidx],'r-',label=lab,lw=0.6)
axs[c].plot(np.arange(sidx-1,70),ens[sidx-1:],'r:',lw=0.6)
axs[c].text(70,ens[-1],bd,va='center',ha='left',color='red',size=8)
axs[0].set_ylabel('Average of $\\langle\\hat{H}\\rangle$')
axs[1].set_xlabel('Sweep')
leg = axs[2].legend(framealpha=1.0)
leg.get_frame().set_linewidth(0.0)
# tit = fig.suptitle('Averages of energy traces for different bond dimensions',weight='heavy')
# top titles
l,b,w,h = axs[0].get_position().bounds
# tit.set_position((tit.get_position()[0],b+h+0.05))
fig.text(l,b+h,' Sweep-dir.:',va='bottom',ha='left',weight='bold')
fig.text(l+w/2,b+h,'Left',va='bottom',ha='center')
# fig.text(l,b+h+0.01,'Sweep sequence',va='bottom',ha='left',size=12,weight='bold')
l,b,w,h = axs[1].get_position().bounds
fig.text(l+w/2,b+h,'Left-Right',va='bottom',ha='center')
l,b,w,h = axs[2].get_position().bounds
fig.text(l+w/2,b+h,'Right',va='bottom',ha='center')
# plt.savefig('...',bbox_inches='tight',dpi=500)
Out[10]:
Text(0.7708333333333335, 0.88, 'Right')
Histogram of results for bond dimension 32, and sweeping leftward¶
In [11]:
measured_points = get_rmse_for_samples(filelocation,filename,np.arange(50),constraint=False,num_sweeps=69,max_bond_dim=32,
sweep_direction='L')
measured_points_constr = get_rmse_for_samples(filelocation,filename,np.arange(50),constraint=True,num_sweeps=69,max_bond_dim=32,
sweep_direction='L')
all_points_sorted = np.sort(measured_points.flatten())
all_points_sorted_constr = np.sort(measured_points_constr.flatten())
print("Average:")
print(f"Without constraint: {np.average(all_points_sorted)}")
print(f"With constraint: {np.average(all_points_sorted_constr)}")
print("Percentage below average")
print(f"Without constraint: {len(all_points_sorted[all_points_sorted<np.average(all_points_sorted)])/len(all_points_sorted)}")
print(f"With constraint: {len(all_points_sorted_constr[all_points_sorted_constr<np.average(all_points_sorted_constr)])/len(all_points_sorted_constr)}")
print("\n50th percentile:")
print(f"Without constraint: {np.percentile(all_points_sorted,50)}")
print(f"With constraint: {np.percentile(all_points_sorted_constr,50)}")
print("\n75th percentile:")
print(f"Without constraint: {np.percentile(all_points_sorted,75)}")
print(f"With constraint: {np.percentile(all_points_sorted_constr,75)}")
Average: Without constraint: 0.0026313240084387688 With constraint: 0.003969939809943597 Percentage below average Without constraint: 0.74 With constraint: 0.764 50th percentile: Without constraint: 0.0008333784619209954 With constraint: 0.0013545194674296171 75th percentile: Without constraint: 0.0028434266393872884 With constraint: 0.0036293316921333445
In [12]:
histrange = (0,0.025)
num_bins = 25*5
fig,axs = plt.subplots(1,2,sharey=True, figsize=(16*cm,8*cm))
fig.subplots_adjust(wspace=0)
ax = axs[0]
ax.grid(lw=0.6)
ax.axvline(np.percentile(all_points_sorted,50),lw=0.6,color="black",ls="solid")
ax.axvline(np.percentile(all_points_sorted,75),lw=0.6,color="black",ls="solid")
ax.axvline(np.average(all_points_sorted),lw=0.6,color="black",ls="dashed")
res = np.histogram(all_points_sorted, bins=num_bins, range=histrange, density=False)
ax.bar(res[1][:-1],res[0],align="edge",width=(histrange[1]-histrange[0])/num_bins,color="red")
ax.set_ylabel("Counts")
ax.set_xlabel("RMSE")
ax.set_xlim(0,ax.get_xlim()[-1])
ax = axs[1]
ax.grid(lw=0.6)
ax.axvline(np.percentile(all_points_sorted_constr,50),lw=0.6,color="black",ls="solid")
ax.axvline(np.percentile(all_points_sorted_constr,75),lw=0.6,color="black",ls="solid")
ax.axvline(np.average(all_points_sorted_constr),lw=0.6,color="black",ls="dashed")
res = np.histogram(all_points_sorted_constr, bins=num_bins, range=histrange, density=False)
ax.bar(res[1][:-1],res[0],align="edge",width=(histrange[1]-histrange[0])/num_bins,color="blue")
ax.set_xlabel("RMSE")
ax.set_xlim(0,ax.get_xlim()[-1])
# plt.savefig("...",bbox_inches="tight")
Out[12]:
(0.0, 0.026250000000000002)
Figure for the sweeping durations¶
In [13]:
# Function for obtaining the averages
def get_sweep_times_averages(filelocation,filename,sample_indices,constraint=None,num_sweeps=None,max_bond_dim=None,
sweep_direction=None,exclude_front=10,exclude_back=10,exclude_top_ratio=0.,exclude_bottom_ratio=0.):
filepath = get_file_from_sample_idx(filelocation,filename,sample_indices[0],
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
num_tries = file['properties'].attrs['num_tries']
if num_sweeps is None:
num_sweeps = file['properties'].attrs['num_sweeps']
elapsed_times = np.empty((len(sample_indices),num_sweeps,num_tries))
for i,sample_idx in enumerate(sample_indices):
filepath = get_file_from_sample_idx(filelocation,filename,sample_idx,
constraint,num_sweeps,max_bond_dim,
sweep_direction,return_full_path=True)
with h5py.File(filepath,'r') as file:
elapsed_times[i,:,:] = np.array(file['data/elapsed_time_sweeps'])
elapsed_times = np.sort(elapsed_times[:,exclude_front:-exclude_back,:].flatten())
num_samples = len(elapsed_times)
elapsed_times = elapsed_times[int(exclude_bottom_ratio*num_samples):(num_samples-int(exclude_top_ratio*num_samples))]
return elapsed_times.mean(),elapsed_times.std()
In [14]:
def legend_title_left(leg):
c = leg.get_children()[0]
title = c.get_children()[0]
hpack = c.get_children()[1]
c._children = [hpack]
hpack._children = [title] + hpack.get_children()
color_red = '#E31033'
color_blue = '#0F33DB'
color_green = '#1BAB36'
color_yellow = '#DBA90F'
num_sweeps = 69
bond_dims = [2,4,8,16,32]
sample_indices = np.arange(50)
sweep_directions = ['L','LR','R']
fig,axs = plt.subplots(1,2,figsize=(16*cm,8*cm),sharey=True)
fig.subplots_adjust(wspace=0)
fitfunc = lambda x,a,b,c,d: a*x**3 + b*x**2 + c*x + d #**2
for ci,constraint in enumerate([False,True]):
ax = axs[ci]
ax.grid()
elapsed_times_per_sweep = np.empty((len(bond_dims),len(sweep_directions),2))
for bdi,bd in enumerate(bond_dims):
for sdi,sd in enumerate(sweep_directions):
elapsed_times_per_sweep[bdi,sdi,:] = get_sweep_times_averages(
filelocation,filename,sample_indices,constraint=constraint,
num_sweeps=num_sweeps,max_bond_dim=bd,
sweep_direction=sd,exclude_front=10,exclude_back=10,exclude_bottom_ratio=0.1,exclude_top_ratio=0.1
)
fit_res_L = curve_fit(fitfunc,bond_dims,elapsed_times_per_sweep[:,0,0])[0]
fit_res_LR = curve_fit(fitfunc,bond_dims,elapsed_times_per_sweep[:,1,0])[0]
fit_res_R = curve_fit(fitfunc,bond_dims,elapsed_times_per_sweep[:,2,0])[0]
print('fit_res_L:',fit_res_L)
print('fit_res_LR:',fit_res_LR)
print('fit_res_R:',fit_res_R)
x = np.arange(len(bond_dims))
xfunc = np.linspace(0,35,100)
ax.set_xticks(bond_dims)#,labels=bond_dims)
ax.plot(xfunc,fitfunc(xfunc,*fit_res_L),'--',color=color_red,alpha=0.5)
ax.errorbar(np.array(bond_dims)-0.35,elapsed_times_per_sweep[:,0,0],yerr=elapsed_times_per_sweep[:,0,1],fmt='.',color=color_red,label='Left')
ax.plot(xfunc,fitfunc(xfunc,*fit_res_LR),'--',color=color_green,alpha=0.5)
ax.errorbar(bond_dims,elapsed_times_per_sweep[:,1,0],yerr=elapsed_times_per_sweep[:,1,1],fmt='.',color=color_green,label='Left-Right')
ax.plot(xfunc,fitfunc(xfunc,*fit_res_R),'--',color=color_blue,alpha=0.5)
ax.errorbar(np.array(bond_dims)+0.35,elapsed_times_per_sweep[:,2,0],yerr=elapsed_times_per_sweep[:,2,1],fmt='.',color=color_blue,label='Right')
# ax.set_yscale('log')
ax.set_xlim(0,34)
ax.set_xlabel('Bond dimension')
axs[0].set_ylim(0,12)
axs[0].set_ylabel('Time/Sweep [s]')
leg = axs[0].legend(title='Sweep dir.:',bbox_to_anchor=(0., 1.0),loc='lower left',ncol=3,
frameon=False,handletextpad=0.05, borderpad=0.,columnspacing=0.5)
legend_title_left(leg)
leg.get_title().set_fontweight('bold')
# plt.savefig('...',bbox_inches='tight',dpi=500)
fit_res_L: [0.00015429 0.00096927 0.0274431 0.11313953] fit_res_LR: [9.70530845e-05 1.13604866e-03 3.17405180e-02 4.50550255e-02] fit_res_R: [ 0.00019596 -0.00043582 0.03848884 0.07430264] fit_res_L: [ 8.12745154e-05 6.48376577e-03 -1.19077065e-02 1.97518201e-01] fit_res_LR: [-2.60272888e-05 8.26066608e-03 -2.26529137e-02 1.64930476e-01] fit_res_R: [1.30259527e-04 4.86120347e-03 7.34782886e-04 1.57778759e-01]
In [ ]:
In [ ]: