Skip to content
Snippets Groups Projects
Commit 603b85c1 authored by Valentin Bruch's avatar Valentin Bruch
Browse files

attempt to created some high-quality plots

parent 7c6fad55
No related branches found
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
final_plots.py 0 → 100644
+ 197
− 0
View file @ 603b85c1
#!/usr/bin/env python3
# Copyright 2022 Valentin Bruch <valentin.bruch@rwth-aachen.de>
# License: MIT
"""
Kondo FRTRG, generate high-quality plots for publication
"""
import scipy.constants as sc
import matplotlib.pyplot as plt
import matplotlib.colors as mplcolors
from matplotlib.widgets import Slider
import argparse
import numpy as np
from scipy.interpolate import bisplrep, bisplev, splrep, BSpline
import settings
from data_management import DataManager, KondoImport
# In this program all energies are given in units of the RTRG Kondo
# temperature Tkrg, which is an integration constant of the E-flow RG
# equations. The more conventional definition of the Kondo temperature is
# G(V=Tk)=G(V=0)/2=e²/h. The ratio Tk/Tkrg is:
TK_VOLTAGE = 3.30743526735
def save_overview(
omega = 16.5372,
vdc_res = 501,
vac_res = 501,
vdc_max = 165.372,
vac_max = 165.372,
method = "mu",
d = 1e9,
xL = 0.5,
solver_tol_rel = 1e-8,
solver_tol_abs = 1e-10,
voltage_branches = 4,
s_g = 1e-4,
s_idc = 5e-6,
s_iac = 5e-6,
filename = "figdata/vdc_vac_omega16.5372_interp.npz",
**kwargs
):
"""
3d plot of ac and dc differential conductance and current as function of Vac and Vdc.
"""
dm = DataManager()
data = dm.list(omega=omega, vdc=None, vac=None, method=method, d=d, xL=xL, solver_tol_abs=solver_tol_abs, solver_tol_rel=solver_tol_rel, voltage_branches=voltage_branches, **kwargs)
data = data.loc[(data.vac < 1.2*vac_max) & (data.vdc < 1.2*vdc_max) & np.isfinite(data.dc_conductance)]
print("Interpolate Gdc", flush=True)
gdc_tck = bisplrep(data.vac, data.vdc, data.dc_conductance, s=s_g, kx=3, ky=3)
print("Interpolate Idc", flush=True)
idc_tck = bisplrep(data.vac, data.vdc, data.dc_current, s=s_idc, kx=3, ky=3)
print("Interpolate Iac", flush=True)
iac_tck = bisplrep(data.vac, data.vdc, data.ac_current_abs, s=s_iac, kx=3, ky=3)
print("done.")
vac_arr = np.linspace(0, min(vac_max, data.vac.max()), vac_res)
vdc_arr = np.linspace(0, min(vdc_max, data.vdc.max()), vdc_res)
gdc_g_interp = bisplev(vac_arr, vdc_arr, gdc_tck).T
idc_interp = bisplev(vac_arr, vdc_arr, idc_tck).T
iac_interp = bisplev(vac_arr, vdc_arr, iac_tck).T
gdc_i_interp = bisplev(vac_arr, vdc_arr, idc_tck, dy=1).T
gac_interp = 2*bisplev(vac_arr, vdc_arr, iac_tck, dx=1).T
np.savez(filename, vac=vac_arr/TK_VOLTAGE, vdc=vdc_arr/TK_VOLTAGE, gdc_g=np.pi*gdc_g_interp, gdc_i=np.pi*gdc_i_interp, gac=np.pi*gac_interp, idc=idc_interp/TK_VOLTAGE, iac=iac_interp/TK_VOLTAGE)
#vac_data, vdc_data = np.meshgrid(vac_arr, vdc_arr)
#with open(filename, 'w') as file:
# file.write("vdc,vac,gdc_g,gdc_i,gac,idc,iac\n")
# np.savetxt(file, np.array([vdc_data/TK_VOLTAGE, vac_data/TK_VOLTAGE, np.pi*gdc_g_interp, np.pi*gdc_i_interp, np.pi*gac_interp, idc_interp/TK_VOLTAGE, iac_interp/TK_VOLTAGE]).reshape((7,-1)).T, fmt="%.9e", delimiter=",")
# #for i in range(vdc_res):
# # file.write('\n')
# # np.savetxt(file, np.array([vdc_data[i]/omega, vac_data[i]/omega, np.pi*gdc_g_interp[i], np.pi*gdc_i_interp[i], np.pi*gac_interp[i], idc_interp[i], iac_interp[i]]).T)
def filter_grid_data(
dm,
omega = 16.5372,
vac_min = 0,
vac_max = 165.372,
vac_num = 101,
vdc_min = 0,
vdc_max = 165.372,
vdc_num = 101,
v_tol = 1e-3,
**kwargs
):
vac_step = (vac_max - vac_min) / (vac_num - 1)
vdc_step = (vdc_max - vdc_min) / (vdc_num - 1)
data = dm.list(omega=omega, vdc=None, vac=None, **kwargs)
data = data.loc[(data.vac >= vac_min - v_tol) & (data.vdc >= vdc_min - v_tol) & (data.vac <= vac_max + v_tol) & (data.vdc <= vdc_max + v_tol) & np.isfinite(data.dc_conductance)]
grid_data = data.loc[(np.abs(((data.vac - vac_min + v_tol) % vac_step) - v_tol) < v_tol) & (np.abs(((data.vdc - vdc_min + v_tol) % vdc_step) - v_tol) < v_tol)]
data = grid_data.copy()
data.sort_values(["vac", "vdc"], inplace=True)
vac_arr = np.linspace(vac_min, vac_max, vac_num)
vdc_arr = np.linspace(vdc_min, vdc_max, vdc_num)
gdc_arr = np.empty((vac_num, vdc_num), dtype=np.float64)
idc_arr = np.empty((vac_num, vdc_num), dtype=np.float64)
iac_arr = np.empty((vac_num, vdc_num), dtype=np.float64)
phase_arr = np.empty((vac_num, vdc_num), dtype=np.float64)
gdc_arr.fill(np.nan)
idc_arr.fill(np.nan)
iac_arr.fill(np.nan)
phase_arr.fill(np.nan)
lower_index = 0
for i, vac in enumerate(vac_arr):
upper_index = data.vac.searchsorted(vac + v_tol)
indices = vdc_arr.searchsorted(data.vdc[lower_index:upper_index] - v_tol)
indices = indices[np.abs(vdc_arr[indices] - data.vdc[lower_index:upper_index]) < v_tol]
gdc_arr[i, indices] = data.dc_conductance[lower_index:upper_index]
idc_arr[i, indices] = data.dc_conductance[lower_index:upper_index]
iac_arr[i, indices] = data.dc_conductance[lower_index:upper_index]
phase_arr[i, indices] = data.dc_conductance[lower_index:upper_index]
lower_index = upper_index
return *np.meshgrid(vdc_arr, vac_arr), gdc_arr, idc_arr, iac_arr, phase_arr
def export_omega5():
omega = 16.5372
dm = DataManager()
# Full overview
vdc, vac, gdc, idc, iac, phase = filter_grid_data(
dm,
omega = omega,
vac_min = 0,
vac_max = 165.372,
vac_num = 101,
vdc_min = 0,
vdc_max = 165.372,
vdc_num = 101,
method = "mu",
d = 1e9,
xL = 0.5,
solver_tol_rel = 1e-8,
solver_tol_abs = 1e-10,
voltage_branches = 4
)
with open("figdata/vdc_vac_omega16.5372.dat", "w") as file:
file.write("vac vdc gdc idc iac")
for i in range(101):
file.write("\n")
np.savetxt(file, np.array([vac[i]/omega, vdc[i]/omega, np.pi*gdc[i], idc[i], iac[i]]).T)
# Lines at constant Vac
for i in range(11):
with open(f"figdata/vdc_vac{i}_omega16.5372.dat", "w") as file:
file.write("vac vdc gdc idc iac\n")
np.savetxt(file, np.array([vac[10*i]/omega, vdc[10*i]/omega, np.pi*gdc[10*i], idc[10*i], iac[10*i]]).T)
# Lines at constant Vdc
for i in range(11):
with open(f"figdata/vac_vdc{i}_omega16.5372.dat", "w") as file:
file.write("vac vdc gdc idc iac\n")
np.savetxt(file, np.array([vac[:,10*i]/omega, vdc[:,10*i]/omega, np.pi*gdc[:,10*i], idc[:,10*i], iac[:,10*i]]).T)
# Zoom to primary Kondo peak
vdc, vac, gdc, idc, iac, phase = filter_grid_data(
dm,
omega = omega,
vac_min = 0,
vac_max = 16.5372,
vac_num = 21,
vdc_min = 0,
vdc_max = 16.5372,
vdc_num = 21,
method="mu",
d=1e9,
xL=0.5,
solver_tol_rel=1e-8,
solver_tol_abs=1e-10,
voltage_branches=4
)
with open("figdata/vdc_vac_omega16.5372_zoom.dat", "w") as file:
file.write("vac vdc gdc idc iac")
for i in range(21):
file.write("\n")
np.savetxt(file, np.array([vac[i]/omega, vdc[i]/omega, np.pi*gdc[i], idc[i], iac[i]]).T)
def prepare_plotly_csv():
dm = DataManager()
# General overview
reduction_dict = dict(omega="omega", vdc="vdc", vac="vac", dc_conductance="g", dc_current="idc", ac_current_abs="iac", ac_current_phase="ac_phase")
data = dm.list(d=1e9, solver_tol_rel=1e-8, solver_tol_abs=1e-10)
data = data.rename(columns=reduction_dict)
data.g *= np.pi
data.omega /= TK_VOLTAGE
data.vdc /= TK_VOLTAGE
data.vac /= TK_VOLTAGE
data.idc /= TK_VOLTAGE
data.iac /= TK_VOLTAGE
data = data.sort_values(["vac","vdc","omega"])
data.to_csv("html/full.csv", columns=reduction_dict.values())
# omega = 5Tk
for omega, name in zip((16.5372, 9.2159791, 5.8206184, 7.1271), ("omega5", "compare_bruhat18a", "compare_bruhat18b", "compare_kogan04")):
reduced = data.loc[np.isclose(data.omega*TK_VOLTAGE, omega) & (data.method == "mu") & (data.voltage_branches == 4)]
reduced.to_csv(f"html/{name}.csv", columns=reduction_dict.values())
# vdc = 0
reduced = data.loc[(data.vdc == 0) & (data.method != "mu")]
reduced.to_csv("html/vdc0.csv", columns=reduction_dict.values())
if __name__ == "__main__":
from sys import argv
globals()[argv[1]]()
plot_pyqtgraph.py
+ 66
− 4
View file @ 603b85c1
......@@ -9,11 +9,18 @@ Kondo FRTRG, generate interactive plots using PyQtGraph
import numpy as np
import argparse
from pyqtgraph.Qt import QtGui
from pyqtgraph import makeQImage
import pyqtgraph.opengl as gl
from matplotlib.pyplot import cm
import pandas as pd
import settings
from data_management import DataManager
#from OpenGL import GL
#GL.glEnable(GL.GL_DEPTH_TEST)
#GL.glShadeModel(GL.GL_SMOOTH)
#GL.glHint(GL.GL_PERSPECTIVE_CORRECTION_HINT, GL.GL_NICEST)
def full_overview(dm, grid=True, size=0.25, xyscale="linear", zscale="log", scale=None, gl_preset="translucent", **parameters):
"""
Show all data for conductance.
......@@ -59,7 +66,7 @@ def full_overview(dm, grid=True, size=0.25, xyscale="linear", zscale="log", scal
w.addItem(sp)
app.exec_()
def fixed_parameter(dm, vac=None, vdc=None, omega=None, scale=80, size=None, grid=False, xyscale="linear", zscale="linear", gl_preset="translucent", **parameters):
def fixed_parameter(dm, vac=None, vdc=None, omega=None, scale=80, size=None, grid=False, xyscale="linear", zscale="linear", gl_preset="translucent", mirror_vdc=True, **parameters):
"""
Show overview of all data where one physical parameter is fixed
"""
......@@ -75,6 +82,11 @@ def fixed_parameter(dm, vac=None, vdc=None, omega=None, scale=80, size=None, gri
assert vac is not None
parameter = "vac"
data = dm.list(vac=vac, vdc=vdc, omega=omega, **parameters)
if mirror_vdc and parameter != "vdc" and xyscale != "log":
data_mirror = data.copy()
data_mirror.vdc *= -1
data = pd.concat((data, data_mirror))
del data_mirror
app = QtGui.QApplication([])
w = gl.GLViewWidget()
w.show()
......@@ -108,6 +120,53 @@ def fixed_parameter(dm, vac=None, vdc=None, omega=None, scale=80, size=None, gri
w.addItem(sp)
app.exec_()
def plot_interpolated(*args, mirror_vdc=True, mirror_vac=False, **parameters):
#data = np.genfromtxt("figdata/vdc_vac_omega16.5372_interp.dat", names=True, delimiter=",")
data = np.load("figdata/vdc_vac_omega16.5372_interp.npz")
app = QtGui.QApplication([])
w = gl.GLViewWidget()
w.show()
#pos = np.array([*np.meshgrid(data["vdc"], data["vac"]), 50*data["gdc_g"]]).reshape((3,-1)).T
#sp = gl.GLScatterPlotItem(pos=pos, size=0.2, color=cm.viridis(data["gdc_g"].flatten()), pxMode=False)
#sp.setGLOptions("additive")
#w.addItem(sp)
#pos_diff = np.array([*np.meshgrid(data["vdc"], data["vac"]), 1000*(data["gdc_i"] - data["gdc_g"])]).T
#sp = gl.GLScatterPlotItem(pos=pos_diff, size=0.2, color=cm.viridis(pos_diff[:,2]/10), pxMode=False)
#sp.setGLOptions("additive")
#w.addItem(sp)
vdc = data["vdc"]
if mirror_vdc:
vdc = np.concatenate((-vdc[:0:-1], vdc))
vac = data["vac"]
if mirror_vac:
vac = np.concatenate((-vac[:0:-1], vac))
gdc = data["gdc_g"]
if mirror_vdc:
gdc = np.concatenate((gdc[:0:-1], gdc), axis=0)
if mirror_vac:
gdc = np.concatenate((gdc[:,:0:-1], gdc), axis=1)
cmap_light = lambda x: cm.viridis((x**.5 - 0.282)/0.718)
cmap_dark = lambda x: 0.1*np.array([[1.,1.,1.,1.]]) + 0.9*cm.viridis((x**.5 - 0.282)/0.718)
sp = gl.GLSurfacePlotItem(x=vdc, y=vac, z=50*gdc, colors=cmap_light(gdc), drawEdges=False)
#sp.translate(-5e-3, -5e-3, -5e-3)
sp.translate(-1e-2, -1e-2, -1e-2)
w.addItem(sp)
for i in range(0, vac.size, 10):
line = gl.GLLinePlotItem(pos=np.array([vdc, vac[i]*np.ones_like(vdc), 50*gdc[:,i]]).T, color=cmap_dark(gdc[:,i]), width=10)
line.setGLOptions("translucent")
w.addItem(line)
for i in range(0, vdc.size, 10):
line = gl.GLLinePlotItem(pos=np.array([vdc[i]*np.ones_like(vac), vac, 50*gdc[i]]).T, color=cmap_dark(gdc[i]), width=10)
line.setGLOptions("translucent")
w.addItem(line)
w.opts["bgcolor"] = (1.,1.,1.,0.)
w.opts["distance"] = 123.5
w.opts["center"] = QtGui.QVector3D(-12.5,0,-10)
arr = w.renderToArray((10000,9700))
makeQImage(arr).save("/tmp/test.png")
app.exec_()
#w.grabFrameBuffer().save('/tmp/test.png')
PRESETS = dict(
default = dict(d=1e9),
vdc0 = dict(vdc=0, d=1e9, scale=30),
......@@ -119,6 +178,7 @@ PRESETS = dict(
# Ω = 5 Tk (Tk defined by G(V=Tk)=e²/h)
omega5 = dict(omega=16.5372, d=1e9),
omega5log = dict(omega=16.5372, xyscale="log", scale=1, d=1e9),
interp = dict(function=plot_interpolated),
)
def main(dm, preset=None, **parameters):
......@@ -127,11 +187,13 @@ def main(dm, preset=None, **parameters):
parameters.update(PRESETS[preset])
except KeyError:
settings.logger.warning("Unknown preset: " + preset)
if "function" in parameters:
parameters["function"]()
return
for name in ("vdc", "vac", "omega"):
if parameters.get(name, None) is not None:
fixed_parameter(dm, **parameters)
break
else:
return
full_overview(dm, **parameters)
def parse():
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment