In [2]:
#from __future__ import division 
%pylab inline
import string
%matplotlib qt 
#%matplotlib notebook
#%matplotlib nbagg
import matplotlib.ticker as ticker
import os

txt = 'run'
#txt = 'OldData/run'
txt2 = '/rho_'
inp ='/input_nrg.dat'


rcParams.update({'font.size': 30})
rc("text", usetex = True)
rc("font", family = "serif")
#rcParams['figure.subplot.hspace'] = 0.3
rcParams['figure.subplot.wspace'] = 0.1
from scipy.interpolate import griddata
from scipy.optimize import curve_fit
from scipy import interpolate
mpl.rcParams['lines.linewidth'] = 2
Gamma = 2.82691*0.01
pi*Gamma


Populating the interactive namespace from numpy and matplotlib
WARNING: pylab import has clobbered these variables: ['griddata']
`%matplotlib` prevents importing * from pylab and numpy
Out[2]:
0.08880999688359521

In [3]:
#read files for folder ini and dot 
#return
def readfile(ini, dot):
        output = txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow_zEQAVG.00.dat'
        
        output = txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow.dat'
        #output = 'NotSoOld_Data/'+txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow.dat'
        #output = 'OldData/'+txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow.dat'
        infile = open(os.path.abspath(output), 'r')
        text = infile.readlines()
        
        #Data managemente
        vec = [];        vec2 = []  
        for x in text:
            #print(list(x))
            a =x.split(' ') ;            vec.append(float(a[-3]))
            vec2.append(float(a[0]))
        #if j%20 == 1 : vec = vec[::-1];
        DOS= array(vec) ; w = array(vec2) 
        return DOS , w
    
zvec = ['0.60','0.80','1.00','1.20','1.40']    
def readfileZ(ini, dot, Ztrick):
        output = txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow_zEQ'+zvec[Ztrick]+'.dat'
        #output = txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow.dat'
        infile = open(os.path.abspath(output), 'r')
        text = infile.readlines()
        #Data managemente
        vec = [];        vec2 = []  
        for x in text:
            #print(list(x))
            a =x.split(' ') ;            vec.append(float(a[-3]))
            vec2.append(float(a[0]))
        #if j%20 == 1 : vec = vec[::-1];
        DOS= array(vec) ; w = array(vec2) 
        return DOS , w    
    
#returns the parameter value of file ini  at the given line    
def readfile_Param(ini, line):
        #output = txt+repr(ini)+txt2 +repr(dot)+'_'+repr(dot)+'_OmegaRhow_zEQAVG.00.dat'
        output = txt+repr(ini)+ inp 
        #output = '/home/cifucito/nrgcode/TwoChNRG/src/Main/Run/Run2DOtM/rho_2_2_OmegaRhow.dat'
        infile = open(os.path.abspath(output), 'r')
        text = infile.readlines()
        return  float(text[line])    
    
    
#Define the log transform. 
#Inputs : w -- vector to transform ,  maxi- max(-log(w)) usually 10
def LogTransform(w, maxi):
    logw = zeros(size(w)) ; logw[w<0] = log10(-w[w<0]) ; logw[w>0] = log10(w[w>0]) 
    logw = -logw ;  #print "b" , logw
    logw[size(logw)/2:] = -logw[size(logw)/2:]+ 2*maxi
    logw =  logw -maxi
    return logw

In [4]:
from mpl_toolkits.axes_grid1 import make_axes_locatable


def LogColorMap( initials , dot , line , f2,ax , Method , maxplot ,PrintAxis , mult):
    W = [] ;T = [] ; Rho_values = [] ; t2=[]
    
    for ini in initials:
        param =  readfile_Param(ini, line) ; t2.append(param)
        
    maxt2 = max(t2) ; mint2 = min(t2); t2 = []
    print maxt2 -mint2
    for ini in initials:
        param =  mult*readfile_Param(ini, line)/(maxt2 - mint2)
        t2.append(param)
        rho , w = readfile(ini, dot) 
        
        #KEY TO READ TAKE THE AVERAGE
        rho1, w1 = readfile(ini+100, dot)
        rho = 0.5*(rho +rho1) 
        
        maxi = max(-log10(w[w>0]))
        logw = LogTransform(w, maxi) ; w = logw
        
        for term in range(size(w)):
            W.append(w[term]) ; T.append(param) ;  Rho_values.append(rho[term])
            
    points = zeros((size(W),2)) ; points[:,0] = array(W) ; points[:,1] = array(T) ; t2 = array(t2); w = array(w)
    Rho_values = array(Rho_values) #/max(Rho_values) 
    grid_x, grid_y = mgrid[min(w):max(w):1000j, min(t2):max(t2):1000j]
    grid_zA = griddata(points, Rho_values, (grid_x, grid_y), method=Method)
    #grid_zA = griddata(points, Rho_values, (grid_x, grid_y), method='linear')
    im = ax.imshow(grid_zA.T, extent=(min(w),max(w),min(t2),max(t2)),vmin=0, vmax=maxplot,cmap="gnuplot2", origin='lower')
    
    divider = make_axes_locatable(ax)
    #cax = divider.append_axes("right", size="5%", pad=0.05);
    #if colorAxis: 
    #ax.set_xlim(-wlims,wlims) ;
    t2 = sort(t2) ; t2 = t2[::2]  
    if PrintAxis :
        paramTicks_scaled = linspace(min(t2),max(t2),3) ; paramTicks = paramTicks_scaled*(maxt2 - mint2)/mult
        yticks(paramTicks_scaled, (str(paramTicks[0]) ,str(paramTicks[1]) ,str(paramTicks[2]) ))
        print 'ParamTicks =' , paramTicks_scaled , paramTicks
        
        ticks = array([-1,-10**(-3) ,-10**(-6) ,10**(-maxi), 10**(-6) ,  10**(-3),1])
        logticks = LogTransform(ticks , maxi)
        xticks(logticks , ('$-1$','$-10^{-3}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-3}$', '$1$'))     
        cbar = f2.colorbar(im)
    #ticks_y = ticker.FuncFormatter(lambda t2, pos: '{0:g}'.format(t2*maxt2/mult))
    #print t2 , ticks_x
    #ax.yaxis.set_major_formatter(ticks_y)
    
    #y1 = mult*0.0025/maxt2
    #ax.axhline(y=y1,linewidth=4, color='r', ls = '--')
    
    #y2 = mult*0.02/maxt2
    #ax.axhline(y=y2,linewidth=4, color='g', ls = '--')
    
    return maxi , im 

def colorMap( initials , dot , line , wlims, f2,ax , Method , maxplot ,PrintAxis , mult):
    W = [] ;T = [] ; Rho_values = [] ; t2=[]
    #mult = 80 ; mult = 0.04 
    for ini in initials:
        param =  readfile_Param(ini, line) ; t2.append(param)
        
    maxt2 = max(absolute(t2)) ; mint2 = min(t2); t2 = []    
    for ini in initials:
        param =  mult*readfile_Param(ini, line)/(maxt2 - mint2)
        t2.append(param)
        rho , w = readfile(ini, dot) 
        
        #KEY TO READ TAKE THE AVERAGE
        rho1, w1 = readfile(ini+100, dot)
        rho = 0.5*(rho +rho1) 
        
        rho = rho[abs(w)<wlims*5];  w = w[abs(w)<wlims*5]/wlims ;
        for term in range(size(w)):
            W.append(w[term]) ; T.append(param) ;  Rho_values.append(rho[term])
    points = zeros((size(W),2)) ; points[:,0] = array(W) ; points[:,1] = array(T) ; t2 = array(t2); w = array(w)
    Rho_values = array(Rho_values) #/max(Rho_values) 
    grid_x, grid_y = mgrid[min(w):max(w):1000j, min(t2):max(t2):1000j]
    grid_zA = griddata(points, Rho_values, (grid_x, grid_y), method=Method)
    #grid_zA = griddata(points, Rho_values, (grid_x, grid_y), method='linear')
    im = ax.imshow(grid_zA.T, extent=(min(w),max(w),min(t2),max(t2)),vmin=0, vmax=maxplot,cmap="gnuplot2", origin='lower')
    
    divider = make_axes_locatable(ax)
    #cax = divider.append_axes("right", size="5%", pad=0.05);
    #if colorAxis: 
    #

    ax.set_xlim(-wlims,wlims) ; t2 = sort(t2) ; t2 = t2[::2]  
    if PrintAxis :
        paramTicks_scaled = linspace(min(t2),max(t2),3) ; paramTicks = paramTicks_scaled*(maxt2 - mint2)/mult
        yticks(paramTicks_scaled, (str(paramTicks[0]) ,str(paramTicks[1]) ,str(paramTicks[2]) ))
        #print 'ParamTicks =' , paramTicks_scaled , paramTicks
        wTicks_scaled = linspace(-1,1,3) ;  wTicks = wTicks_scaled*wlims
        xticks(wTicks_scaled, (str(wTicks[0]) ,str(wTicks[1]) ,str(wTicks[2])))                           
        #print 'wTicks =' , wTicks_scaled , wTicks    
        
    #ticks_y = ticker.FuncFormatter(lambda t2, pos: '{0:g}'.format(t2*maxt2/mult))
    #print t2 , ticks_x
    #ax.yaxis.set_major_formatter(ticks_y)
    
    #y1 = mult*0.0025/maxt2
    #ax.axhline(y=y1,linewidth=4, color='r', ls = '--')
    
    #y2 = mult*0.02/maxt2
    #ax.axhline(y=y2,linewidth=4, color='g', ls = '--')
    #return max(Rho_values) , im ]
    return im 


#Plots 2 or 4 Dos plots
#Input : LogPlot (boolean) : true - logPlot , false - Normal plot
#        numplots (1,2) : 1 - plots only dots 0 and 1 , 2 - plots 4 dots
#        initials - array of file numbers
#        line - line of changing parameter 
#        wlims - defines the w limits if logPlot = false 
#        paramName - name of changing param 
#        Method: Method of approximation - 'Linear ' , 'Nearest'
#        MaxRho: Maximum Permited DOS 
#        mult : Allows to change the dimensions of the plot if required
def color4Dots( LogPlot, numPlots,  initials ,  line , wlims ,paramName, Method , MaxRho , mult ):  
    f2, axarr = subplots(numPlots, 2, sharex='col', sharey='row', figsize=(16, 5*numPlots))
    maxplot = 0 ; dots = range(2,4) ;  dots = range(2*numPlots) 
    t2 = [] 
    for ini in initials:
        param =  readfile_Param(ini, line) ;  t2.append(param)
    maxt2 = (max(t2) -min(t2)) ; t2 = sort(mult*array(t2)/maxt2 ) ;  t2 = t2[::2] 
    print maxt2 , t2
    #for dot in dots:
    #    for ini in initials:
    #        rho , w = readfile(ini, dot) 
    #        maxrho = max(rho[abs(w)<wlims]) 
    #        #print maxrho , ini , dot
    #        maxplot = maxrho if maxrho>maxplot else maxplot 
    for dot in dots:
        spin = '\uparrow'
        if dot%2 == 1 : spin = '\downarrow'; 
        ax = axarr[dot-dots[0]] if numPlots == 1 else axarr[dot/2 , dot%2 ]
        
        if LogPlot:
            maxi , im = LogColorMap( initials , dot , line , f2, ax, Method, MaxRho , False, mult)
            ticks = array([-1,-10**(-3) ,-10**(-6) ,10**(-maxi), 10**(-6) ,  10**(-3),1])
            wTicks_scaled = LogTransform(ticks , maxi)
            wTicksLabels = ['$-1$','$-10^{-3}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-3}$', '$1$']
        else:
            im = colorMap(initials , dot , line , wlims,f2, ax, Method ,MaxRho  , False , mult)
            wTicks_scaled = linspace(-1,1,3) ;        wTicks = wTicks_scaled*wlims
            wTicksLabels = ['$'+str(wTicks[0])+'$', '$'+str(wTicks[1])+'$', '$'+str(wTicks[2])+'$']
            
        paramTicks_scaled = linspace(min(t2),max(t2),3) ;
        paramTicks = paramTicks_scaled*maxt2/mult    
        print  paramTicks_scaled, paramTicks
        
        setp(axarr, xticks = wTicks_scaled, xticklabels= wTicksLabels ,
        yticks=paramTicks_scaled ,yticklabels=['$'+str(paramTicks[0])+'$', '$'+str(paramTicks[1])+'$', '$'+str(paramTicks[2])+'$'] )
        
        #setting ticks
        
        title = '$\\rho_{1,2'+ spin +'}$' if numPlots == 1 else '$\\rho_{'+str(dot/2 + 1)+ spin +'}$'
        ax.set_title(title)
        if dot >= (numPlots*2-2) : ax.set_xlabel('$\omega/D'+ '$')
        if dot%2 == 0 : ax.set_ylabel('$'+paramName+'$') 
            
    cbar = f2.colorbar(im, ax=axarr.ravel().tolist() )
    cbar.set_ticks([0, 0.5*MaxRho, MaxRho])
    cbar.set_ticklabels(['$0$', '$'+str(0.5*MaxRho)+'$', '$>'+str(MaxRho)+'$'])
    #cbar.ax.set_yticklabels(['< -1', '0', '> 1'])
    #f2.subplots_adjust(right=0.8) ticks=[0, 0.5*MaxRho, MaxRho]
    #cbar_ax = f2.add_axes([0.85, 0.15, 0.05, 0.7])
    #f2.colorbar(im, cax=cbar_ax)

In [65]:


In [131]:
dot=4
line = 13
Method = 'linear' #; Method = 'nearest'
initials = array([450,451,452,453,454,455,456,457,458,459])
#initials = range(400,409)

#initials = range(100,109)
f2, ax = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 8))

maxplot = 0
for ini in initials:
     rho , w = readfile(ini, dot) ;     maxrho = max(rho[abs(w)<wlims*2]) ;    maxplot = maxrho if maxrho>maxplot else maxplot 
LogColorMap( initials , dot , line , f2, ax, Method, maxrho , True , 10)
print maxrho  , maxplot


0.0225
ParamTicks = [ 0.          4.44444444  8.88888889] [ 0.    0.01  0.02]
17.8445308803 528.303036741

In [16]:
wlims = 0.1
dot=4
line = 18
Method = 'linear' #; Method = 'nearest'
initials = array([450,451,452,453,454,455,456,457,458,459])

initials = range(100,109)
f2, ax = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 8))

maxplot = 0
for ini in initials:
     rho , w = readfile(ini, dot) ;     maxrho = max(rho[abs(w)<wlims*2]) ;    maxplot = maxrho if maxrho>maxplot else maxplot 
colorMap( initials , dot , line , wlims,f2, ax, Method, maxrho , True , 1.5)
print maxrho  , maxplot

#color4Dots( initials ,  line , wlims)


ParamTicks = [ 0.    0.75  1.5 ] [ 0.   0.1  0.2]
wTicks = [-1.  0.  1.] [-0.1  0.   0.1]
1.03022887428 35.7754363602

In [159]:
wlims = 0.05
Method = 'linear' #; Method = 'nearest';  
line = 19
xtitle = [None]*6
xtitle[0] = 't_1/D = t_2/D' ; xtitle[1] = '\epsilon_{d2}/D'; xtitle[2] = '\epsilon_{d2}/D = \\frac{-U_2}{2D}';xtitle[3] = 't_2/D'
xtitle[5] = 't_{dots}/D ' 
rcParams.update({'font.size': 30})
initials = array([450,451,452,453,454,455,456,457,458,459])
initials = array([460,461,462,463,464,465,466,467,468,469])
#initials = array([470,471,472,473,474,475,476,477,478,479])
#initials = array([430,431,432,433,434])
initials = range(400,409)
#initials = array([420,421,422,423,424,425,426,427,428,429])
initials = range(430,439)
#initials = range(410,415)

initials = array([100,130,131,132,133,134])
initials = range(140,144)
#initials = range(110,119)
#initials = range(120,129)
#initials = range(420,429)
#initials = range(430,434)


#Normal 2D
#color4Dots(False, 2,initials ,  line , wlims, xtitle[5], Method ,3,2)

#Log 2D
color4Dots(True, 2,initials ,  line , wlims, xtitle[5], Method ,1,12)


0.9 [  4.18878933  12.18878933]
0.9
[  4.18878933   8.18878933  12.18878933] [ 0.3141592  0.6141592  0.9141592]
0.9
[  4.18878933   8.18878933  12.18878933] [ 0.3141592  0.6141592  0.9141592]
0.9
[  4.18878933   8.18878933  12.18878933] [ 0.3141592  0.6141592  0.9141592]
0.9
[  4.18878933   8.18878933  12.18878933] [ 0.3141592  0.6141592  0.9141592]

In [19]:


In [56]:


In [ ]:


In [54]:
def plotFermi(initials,line,limsx, limsy, Param_name):
    f2, axarr2 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 8))
    #plot the relation between both DOS (up/Down)
    f3, axarr3 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 8))
    t2 = []    
    for ini in initials: 
        t2.append(readfile_Param(ini, line))
    t2 = array(t2)
    
    for i in range(2):
        vec0 = [] ;         maxi = []
        #for j in range(40,46):
        for j in initials:
            rho , vec2 = readfile(j, i)
            
            #KEY TO READ TAKE THE AVERAGE
            rho1, w1 = readfile(j+100, i)
            vec = 0.5*(rho +rho1) 
            
            vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
            vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]
            vec0.append(0.5*(vec0mas + vec0menos)) ; maxi.append(max(vec))
        spin = '\uparrow' if i%2 == 0 else '\downarrow';  
        title = '$\\rho_{'+str(i/2 +1 )+ spin +'}$'    
        axarr2.plot(t2,vec0 ,'o', label= title + '$(0)$')
        if i%2 == 0 : 
            new0 = vec0
        else: 
            ydata = array(new0)/array(vec0) ;  t2_new = sort(t2) ;  ydata_new = ydata[argsort(t2)] ; 
            axarr3.plot(t2_new,ydata_new, label= 'Dot' + str(i/2 + 1),  marker="o" , markersize=10, linewidth=3.0)
            ydata_new = ydata[argsort(t2)] ; 
            #axarr3.plot(t2_new,ydata_new,'-', label= 'Dot' + str(i/2 + 1) )
            #tnew =  np.linspace(t2[0], t2[-1], 50) ; f = interpolate.interp1d(t2, ydata, kind='linear')
            #axarr3.plot(tnew,f(tnew),'-' )

    axarr2.set_xlabel('$'+Param_name+'$')
    axarr2.set_ylabel('$\\rho(0)$')
    axarr2.legend(fontsize=20 , loc='upper right' )
    #axarr2.set_xlim([-0.25,-0.05])
    #axarr2.set_ylim([0.0,4.5])
    #axarr3.set_xlabel('$\epsilon_{d2}$',fontsize=25)
    axarr3.set_xlabel('$'+Param_name+'$')
    axarr3.set_ylabel('$\\rho_{\uparrow}(0)/\\rho_{\downarrow}(0)$')     
    axarr3.legend(fontsize=27, loc='upper right')  
    axarr3.set_xlim([t2[0],t2[-1]])
    t2 = sort(t2) 
    paramTicks = linspace(t2[1],t2[-2],3) 
    axarr3.set_xticks(paramTicks) ;    axarr3.set_xticklabels(['$'+str(paramTicks[0])+'$', '$'+str(paramTicks[1])+'$', '$'+str(paramTicks[2])+'$'])
    axarr3.set_yticks(array([1,2,3,4,5])) ;  axarr3.set_yticklabels( ['$1.0$', '$2.0$', '$3.0$', '$4.0$' ,'$5.0$'])
    #axarr3.set_yticks(linspace(1,1.5,2))
    #axarr3.set_yticklabels(['$1$', '$1.5$', '$2$'])
    #axarr3.set_xticks(linspace(-2,2,5), ('$-2$','$-1$', '$0$','$1$', '$2$'))
    #axarr3.set_ylim([0,3])

In [55]:
xtitle = 't_1/D = t_2/D' ; xtitle = '\epsilon_{d2}/D'; #xtitle = '\epsilon_{d2}/D = \\frac{-U_2}{2D}';# xtitle = 't_2/D'
rcParams.update({'font.size': 45})
initials = array([450,451,452,453,454,455,456,457,458,459])
initials = array([460,461,462,463,464,465,466,467,468])
initials = range(400,409)
initials = range(420,429)
initials = range(110,119)
initials = range(120,129)
#initials = array([420,421,422,423,424,425,426,427,428,429])
#initials = array([470,471,472,473,474,475,476,477,478,479])
plotFermi(initials,4,0.1,4,xtitle)

In [ ]:


In [28]:
#Define the log transform. 
#Inputs : w -- vector to transform ,  maxi- max(-log(w)) usually 10
def LogTransform(w, maxi):
    logw = zeros(size(w)) ; logw[w<0] = log10(-w[w<0]) ; logw[w>0] = log10(w[w>0]) 
    logw = -logw ;  #print "b" , logw
    logw[size(logw)/2:] = -logw[size(logw)/2:]+ 2*maxi
    logw =  logw -maxi
    return logw 
        #    return logVector
def LogPlotDOS(initials , dot , param , param_name):
    for i in initials:
        rho, w = readfile(i, dot)
        #rho1, w1 = readfile(i+100, dot)
        #rho = 0.5*(rho +rho1) ; w = 0.5*(w+w1)
        #maxi = max(-log(w[w>0])) ; print maxi
        maxi = max(-log10(w[w>0]))
        logw = LogTransform(w , maxi)
        ticks = array([-1,-10**(-3) ,-10**(-6) ,10**(-maxi), 10**(-6) ,  10**(-3),1])
        logticks = LogTransform(ticks , maxi)
        xticks(logticks , ('$-1$','$-10^{-3}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-3}$', '$1$'))
        #rhoDw, w3 = readfile(i, dot+1)
        #rhoDw1, w4 = readfile(i+100, dot+1)
        #rhoDw = 0.5*(rhoDw +rhoDw1) 
        #logw = zeros(size(w)) ; logw[w<0] = log10(-w[w<0]) ; logw[w>0] = log10(w[w>0]) 
        #logw = -logw ; 
         #print "b" , logw
        spin = '\uparrow' if dot%2==0 else '\downarrow'
        pr = readfile_Param(i, param)
        lb = '$'+spin +','+ param_name +'= '+str(pr)+ 'D$ ' 
        style = '-' if dot%2==0 else '--'
        #if pr == 0:
        #    lb = '$'+param_name +'= '+ str(pr) + 'D$ ' 
        #    plot(logw, rho,label = lb ,color = 'k', marker = 'o',  ls = 'None' , markersize=5 )
        #else:    
        plot(logw, rho, label = lb , ls = style)
        #    mk = '^' if dot%2==0 else 'v'
        #    #mk = "*" if dot%2==0 else 'H'
        #    if i%2==0:
        #        cl = 'r'# if dot%2==0 else 'y'
        #    else:    
        #        cl = 'b'# if dot%2==0 else 'g'
        #plot(logw, rho,label = lb ,color = cl, marker = mk,  ls = 'None' , markersize=-20/log10(pr))
        #plot(logw, rho, label = lb , ls = style)
        #plot(logw, log2(rho/rhoDw), label = lb , ls = style)
        #plot(logw, rhoDw, label = lb , ls = style)
        #plot(logw, rhoDw, label = lb + 'dw' , ls = style)
            #logw[size(logw)/2:]= logw[size(logw)/2:]+1
    return logticks
    
#print logticks , str(readfile_Param(i, param))
xlabel('$\omega/D$')
ylabel('$\\rho_{1,2}/D$')
#title('Dot 1  $(t_1=t_2=0.005)$')

initials = [458 , 451 , 450] ; 
#initials = [451,458 ] ;
#initials = [460,466,468 ] ;
#initials = [430,434 ] ;
#initials = [400,403 ,407 ] ;
initials = [410 , 411 , 414] ; 
initials = [32,33]
#initials = [22]
param = 21
xtitle = 't_{1,2}/D' ; #xtitle = '\epsilon_{d2}/D'#; xtitle = '\epsilon_{d2}/D = \\frac{-U_2}{2D}';# xtitle = 't_2/D'
LogPlotDOS(initials, 0,param , xtitle)
LogPlotDOS(initials, 1,param , xtitle)
#logticks=LogPlotDOS(initials[:-1] , 1,param , xtitle)
#LogPlotDOS(551 , param)
legend(fontsize=18,loc='lower right')
#title('Dot 2 ($t_1=t_2=0.0025D$)')


Out[28]:
<matplotlib.legend.Legend at 0x7f557bb1ff90>

In [10]:
xticks(logticks , ('$-10^{-2}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-2}$'))
legend(fontsize=16)


---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-10-65f37a917c47> in <module>()
----> 1 xticks(logticks , ('$-10^{-2}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-2}$'))
      2 legend(fontsize=16, loc='lower right')

NameError: name 'logticks' is not defined

In [ ]:


In [ ]:


In [23]:
def plotDOS(initials,line,Param_name,limsx, limsy):

    #Param_name = '\epsilon_2 = \\frac{-U_2}{2}'
    #Param_name = 't_2'
    
    f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
    #plot the DOS at the FE
    t2 = []
    for ini in initials: 
        t2.append(readfile_Param(ini, line))
    t2 = array(t2) 
    
    #suptitle('Majorana Connected to the First QD')
    #print initials
    
    for i in range(4):
        vec0 = []
        maxi = []
        #for j in range(40,46):
        for j in initials:
            param = t2[where(initials == j)] ; #param =  param[0]
            
            #Read file
            vec , vec2 = readfile(j, i)
            #plot DOS
            axarr[i/2 , i%2 ].plot(vec2,vec, '-',label= '$' + Param_name + '= '+str(param) + '$') 

        #Set graphic font
        if i%2 == 1: axarr[i/2 , i%2 ].legend(fontsize=16)
        #axarr[i/2 , i%2 ].set_ylim([0,max(maxi)])    
        #axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
        spin = '\uparrow' 
        if i%2 == 1 : spin = '\downarrow';  
        title = '$\\rho_{'+str(i/2 +1 )+ spin +'}$'
        axarr[i/2 , i%2 ].set_title(title)    
        axarr[i/2 , i%2 ].set_xlim([-limsx,limsx])
        axarr[i/2 , i%2 ].set_ylim([0,limsy])
        
        if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\epsilon $',fontsize=26)
        if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('$\\rho_{'+str(i/2 +1 )+ '}$',fontsize=26)
xtitle = [None]*23           
xtitle[12] = 't_1/D = t_2/D' ; xtitle[16] = '\epsilon_{d2}/D'; xtitle[15] = '\epsilon_{d2}/D = \\frac{-U_2}{2D}';xtitle[13] = 't_2/D'

xtitle[17] = 't_{dots}/D ' ; xtitle[18] = '\epsilon_m ' ;xtitle[19] = '\phi_1' ;  xtitle[20] = '\phi_2' ;xtitle[21] = '\eta_1' ;
xtitle[22] = '\eta_1'      

line = 21
plotDOS([22],line,xtitle[line],1, 1)

In [ ]:


In [ ]:


In [ ]:


In [260]:
from scipy.signal import argrelextrema
i = 451 ; dot = 1
def findLocalMax(i, dot):
    rho, w = readfile(i, dot) ;rho1, w1 = readfile(i+100, dot)
    rho = 0.5*(rho +rho1) ; w = 0.5*(w+w1)

    arg = array(argrelextrema(rho[1::2], greater))
    w1 = -array(w[arg])
    a = w1>0.0001 
    b = w1<0.5
    w1= w[arg[a&b]] ; rho1 = rho1[arg[a&b]]
    return w1,rho1

In [242]:


In [ ]:


In [190]:
xtitle = 't_1/D = t_2/D' ; xtitle = '\epsilon_{d2}/D';# xtitle = '\epsilon_{d2}/D = \\frac{-U_2}{2D}'; xtitle = 't_2/D'
initials = array([56,57,58,59,60,256,257,258,259,260])
initials = array([261,262,263,264,265,266,267,268,269])
initials = array([8,10])


rho, w = readfile(8, 0)
maxi = max(rho)
rho2, w2 = readfile(8, 1)
Gamma = 2.82691*0.01
print pi , Gamma
plot( w/Gamma , rho/maxi ,  label = 'QSz $\\rho_\uparrow$' ,color='b', marker="D" , linewidth=2.0)
plot(w2/Gamma , rho2/maxi , label = '$QSz \\rho_\downarrow$',color='g', marker="D" ,linewidth=2.0)

rho, w = readfile(10, 0)
maxi = max(rho)
rho2, w2 = readfile(10, 1)
Gamma = 2.82691*0.01
print pi , Gamma
plot( w/Gamma , rho/maxi ,  label = 'NupPdn $\\rho_\uparrow$' ,color='b', marker="s" , linewidth=2.0)
plot(w2/Gamma , rho2/maxi , label = 'NupPdn $\\rho_\downarrow$',color='g', marker="s" ,linewidth=2.0)

legend(fontsize=26)
xlim(-2,2)
xlabel('$\epsilon/\Gamma$') ; ylabel('Density of States $(\pi \Gamma \\rho)$')
yticks(linspace(0,1,3), ('$0$', '$0.5$', '$1$'))
xticks(linspace(-2,2,5), ('$-2$','$-1$', '$0$','$1$', '$2$'))
#plotDOS(initials,16,1,4)
print maxi


3.14159265359 0.0282691
3.14159265359 0.0282691
2.37161046419

In [55]:
def plotPH_Sym(initials , param , mini,maxi):
    f2, ax = subplots(size(initials), 1, sharex='col', sharey='row', figsize=(10, 4*(size(initials))))
    #suptitle('Majorana Connected to the First QD')
    for i in initials:
        num = where(initials == i)[0] ; num = num[0]
        for j in array([0,2]):
            dosUp , wUp = readfile(i, j)
            dosDown , wDown = readfile(i, j+1)
            #plot(wUp , dosUp, wUp , dosDown)
            vec1 = dosUp[ wUp > 0] + dosDown[ wDown > 0]
            vec2 = dosUp[ wUp < 0] + dosDown[ wDown < 0]
            #vec1 = dosUp[ wUp > 0] - dosUp[ wUp < 0][::-1] 
            #vec2 = dosDown[ wDown > 0][::-1] - dosDown[ wDown < 0]

            #print vec1 - vec2[::-1]
            #plot(wUp[ wUp > 0],vec1, wUp[ wUp < 0],vec2 )
            #ax.semilogx( wUp[ wUp > 0] , vec1/vec2[::-1] , 'o',label = 'D '+ str(j/2+1)+', i '+str(i))
            #if num%3 == 0 :
            #ax[num].semilogx( wUp[ wUp > 0] , vec2[::-1]/vec1, 'o',label = 'Dot '+str(j/2+1))
            #else :
            ax[num].loglog( wUp[ wUp > 0] , vec1/vec2[::-1], 'o',label = 'Dot '+str(j/2+1))
            #plot( range(size(vec1))[::-1] , vec1/vec2[::-1] , 'o',label = 'Dot '+ str(j/2 +1))
            #print vec1[-1], vec2[0], vec1[-1]/vec2[0] 
            #plot( log(wUp[ wUp > 0]), vec1 , 'o',label = ''+ str(j/2))
            #plot(  (-1)*log(wUp[ wUp > 0]) , vec2[::-1] , 'o',label = ''+ str(j/2))
            #xlim([-10,10])
            #ylim([1.0-0.00002,1.0+0.00002])
            #ax[num].set_title('$\Gamma =' + str(readfile_Param(i, 3))+', N_s = ' +repr(readfile_Param(i, 1))+', \epsilon_{d2} = ' +repr(readfile_Param(i, 16))+'$'  , fontsize=22)
            ax[num].set_title(' h_z = ' +repr(readfile_Param(i, param))+'$'  , fontsize=22)
        
        #ax[num].text(1, 1, r'$E=mc^2$', fontsize=15) 
        #ax[num].set_ylim(mini,maxi)
        ax[num].set_ylabel('$\\frac{\\rho_\uparrow(\omega^+)+ \\rho_\downarrow(\omega^+)}{\\rho_\uparrow(\omega^-) +\\rho_\downarrow(\omega^-)}$')
            #print vec1 - vec2[::-1]
    ax[num].legend(fontsize=20 , loc='upper left' )    
    ax[num].set_xlabel( '$\omega$' )  
                

        
def plotDOS(initials,line,limsx, limsy):
    Param_name = 't_1 = t_2'
    Param_name = '\epsilon_2'
    #Param_name = '\epsilon_2 = \\frac{-U_2}{2}'
    #Param_name = 't_2'
    
    f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
    #plot the DOS at the FE
    t2 = []
    for ini in initials: 
        t2.append(readfile_Param(ini, line))
    t2 = array(t2) 
    
    #suptitle('Majorana Connected to the First QD')
    #print initials
    
    for i in range(4):
        vec0 = []
        maxi = []
        #for j in range(40,46):
        for j in initials:
            param = t2[where(initials == j)] ; param =  param[0]
            
            #Read file
            vec , vec2 = readfile(j, i)
            #plot DOS
            axarr[i/2 , i%2 ].plot(vec2,vec, '-',label= '$' + Param_name + '= '+str(param) + '$') 

        #Set graphic font
        if i%2 == 1: axarr[i/2 , i%2 ].legend(fontsize=16)
        #axarr[i/2 , i%2 ].set_ylim([0,max(maxi)])    
        #axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
        spin = '\uparrow' 
        if i%2 == 1 : spin = '\downarrow';  
        title = '$\\rho_{'+str(i/2 +1 )+ spin +'}$'
        axarr[i/2 , i%2 ].set_title(title)    
        axarr[i/2 , i%2 ].set_xlim([-limsx,limsx])
        axarr[i/2 , i%2 ].set_ylim([0,limsy])
        
        if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\epsilon $',fontsize=26)
        if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('$\\rho_{'+str(i/2 +1 )+ '}$',fontsize=26)




    
def logplotDOS(initials , line):
    #plot the 4 DOS
    f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
    #plot the DOS at the FE
    f2, axarr2 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
    #plot the relation between both DOS (up/Down)
    f3, axarr3 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
    
    t2 = []
    for j in initials: 
        t2.append(readfile_Param(j, line))
    t2 = array(t2) 
    for i in range(4):
        vec0 = []
        #for j in range(40,46):
        for j in initials:
            param = t2[where(initials == j)] ; param =  param[0]
            #Read file
            vec , vec2 = readfile(j, i)
            #get 0 value of the DOS
            if j%10 ==0 : vec = vec[::-1];
            vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
            vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]

            vec0.append(0.5*(vec0mas + vec0menos))

            #plot DOS
            #axarr[i/2 , i%2 ].plot(vec2,vec, label= '$\epsilon_{d2} ='+str(-0.25+(j-initials[0])*0.05 )+'$') 
            axarr[i/2 , i%2 ].semilogx(vec2[ vec2 > 0],vec[ vec2 > 0],'o',  label= '$'+ str()+'+$') 
            axarr[i/2 , i%2 ].semilogx(-vec2[ vec2 < 0],vec[ vec2 < 0],'o' ,label= '$'+ str(param)+',-$') 
        #Set graphic font
        if i%2 == 1: axarr[i/2 , i%2 ].legend(fontsize=12)
        #axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
        spin = '\uparrow' 
        if i%2 == 1 : spin = '\downarrow';  
        title = '$\\rho_{'+str(i/2 +1 )+ spin +'}$'
        print title
        axarr[i/2 , i%2 ].set_title(title)    
        #axarr[i/2 , i%2 ].set_ylim([-0.0,9.0])
        if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\epsilon $',fontsize=25)
        if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('$\\rho_{'+str(i/2 +1 )+ '}$',fontsize=25)

        #t2 = -0.25 + double(initials-initials[0]) * 0.05  
        #t2[-1] = -0.06 ; t2[-2] = -0.08 
        axarr2.plot(t2,vec0 ,'o', label= title + '$(0)$')


        if i%2 == 0 : 
            new0 = vec0
        else: 
            print array(vec0)/array(new0) 
            axarr3.plot(t2,array(new0)/array(vec0) ,'o', label= 'Dot' + str(i/2 + 1))

    axarr2.set_xlabel('$\epsilon_{d2}$',fontsize=25)
    axarr2.set_ylabel('$\\rho(0)$',fontsize=25)
    axarr2.legend(fontsize=24 , loc='upper right' )
    #axarr2.set_xlim([-0.25,-0.05])
    #axarr2.set_ylim([0.0,4.5])

    axarr3.set_xlabel('$\epsilon_{d2}$',fontsize=25)
    axarr3.set_ylabel('$\log_2(\\rho_\uparrow(0)/\\rho_\downarrow(0))$',fontsize=25)        
    axarr3.legend(fontsize=20, loc='upper right')        
    #axarr3.set_xlim([-0.25,-0.05])
    #axarr3.set_ylim([0.0,1.8])
    #save = directory +'/'+txt2 +repr(i)+'_'+repr(i)+'_OmegaRhow.png'
    #savefig(save)    

zvec2 = ['0.6','0.8','1','1.2','1.4']    
def readCost(ini, Ztrick):
        output = txt+repr(ini)+'/output_DMNRG_zEQ'+zvec2[Ztrick]+'_Dir'+repr(ini)+'.txt'
        #output = txt+repr(ini)+'/output_NRG_Code_DM_NRG_Dir'+repr(ini)+'.txt'
        infile = open(os.path.abspath(output), 'r')
        text = infile.readlines()
        w= []
        costi = []
        for x in text:
                #print(list(x))
            a=x.split('(')
            if a[0]=='RhoCosti0_0':
                omega = a[1].split('=')[1] ; omega = omega.split(' ')[1] 
                cost = a[1].split('=')[2] 
                w.append(float(omega)) ; costi.append(float(cost))
                #print float(omega) ,cost
                #Data managemente
        w = array(w) ; costi = array(costi)  
        costi = costi[argsort(w)] ;         w = w[argsort(w)]
        return w , costi

In [17]:
initials = array([60])
#initials = array([70,71,72,73,74])
#initials = initials[::2]
initials = array([109,119,124])
initials = array([105 ,106 ,107])
initials = array([160,161,162,163,164])
initials = array([160,165])
initials = array([61,62,63,64,65])
initials = array([231,232,233,234,235,236,237,238,239])
#initials = array([70,71,72,73,74])
#initials = array([56,57,58,59,60])
initials = array([451,551])

#initials = array([56,57,58,59,60,256,257,258,259,260])
#initials = array([35,36,37,38,39])
#initials = array([75,76,77,78,79])
#initials = array([80,81,82,83,84])

#initials = initials[::2]
#plotDOS(initials,line,limsx, limsy)
#initials - vector with initial data  files , line - changing parameter , limsx,limsy - plot limit axis 
rcParams['figure.subplot.hspace'] = 0.15
rcParams['figure.subplot.wspace'] = 0.0
plotDOS(initials,13,0.1,4)

In [58]:
initials = array([70,71,72,73,74])
#initials = array([70])
#initials = array([6])
initials = array([2002])
#print readfile_Param(1001, 1)
#print readfile_Param(1001, 3)
#initials = array([2011,2012,2013])
initials = array([2012,2013,2014])
initials = array([2020 , 2021, 20])
#initials = array([70,71,72,73,74])
initials = array([56,57,58,59,60])
initials = array([451])
initials = array([451,551])
#initials = initials[::2]

plotPH_Sym(initials,11,-10,10)
suptitle('Majorana - DQD $(t_1=t_2 =0.005)$')


Out[58]:
<matplotlib.text.Text at 0x7ff22256afd0>

In [31]:
readfile_Param(2013, 3)


Out[31]:
0.088

In [14]:
initials = array([85,86])
#initials = initials[::2]
logplotDOS(initials,17)


$\rho_{1\uparrow}$
$\rho_{1\downarrow}$
[ 0.75854088  0.54638411]
$\rho_{2\uparrow}$
$\rho_{2\downarrow}$
[ 0.51861326  0.50936563]

In [153]:
#param = t2[where(initials == j)] ; param =  param[0]

ini =145
dot = 0
rcParams['figure.subplot.wspace'] = 0.2
f, axarr = subplots(1, 2, figsize=(14, 6))
params = [21,23,25,27,29,31]
#for ini in range(size(initials)):
for z in range(size(zvec)):
    vec , vec2 = readfileZ(ini, dot,z)
    vec1 , vec21 = readfileZ(ini, dot+1,z)
    #vec , vec2 = readfile(initials[ini], dot)
    #axarr[0].semilogx(vec2[vec2>0],vec[vec2>0]/vec[vec2<0][::-1])
    #plot(vec2,vec,'o',label= '$Z ='+ zvec[z] + '$')
    axarr[0].semilogx(vec2,vec,'-',label= '$Z ='+ zvec[z] + '$')
    axarr[1].plot(vec21,vec1,'-',label= '$Z ='+ zvec[z] + '$')
    #axarr[1].plot(vec2,vec, label= '$N ='+ repr(params[ini]) + '$')
    #xlabel('$\omega$',fontsize=25)        
    #axarr[1].set_xlabel('$\omega$',fontsize=25)
    #ylabel('$\\rho$',fontsize=25)        
    #axarr[1].set_ylabel('$\\rho(\omega)$',fontsize=25)   
    legend(fontsize=18 , loc='upper right' )
    #axarr[0].set_ylim(0,2) 
    axarr[1].set_xlim(-1 , 1) 
    axarr[0].set_xlim(-1 , 1)

In [ ]:


In [45]:
rcParams['figure.subplot.wspace'] = 0.2
f, axarr = subplots(1, 2, figsize=(14, 6))
#params = [0.001 , 0.01,3]
params = zvec2 
print params
#print readCost(ini, 0)
#for ini in range(size(initials)):
ini =26
for z in range(size(zvec2)):
        w , costi = readCost(ini, z)
        #w , costi = readCost(initials[ini], z)
        axarr[0].semilogx(w[w>0],costi[w>0]/costi[w<0][::-1])
        axarr[1].plot(w,costi,label= 'Z =' + zvec2[z])
        axarr[1].plot(w,costi,label= 'Z =' + zvec2[z])
        #axarr[1].plot(w,costi,'o',label= '$N ='+ repr(params[ini]) + '$')    
        axarr[1].plot(w,costi,'-') 
        axarr[0].set_xlabel('$\omega$',fontsize=25)        
        axarr[1].set_xlabel('$\omega$',fontsize=25)
        axarr[0].set_ylabel('$cost(\omega +)/cost(\omega -)$',fontsize=25)        
        axarr[1].set_ylabel('$cost(\omega)$',fontsize=25)   
        axarr[1].legend(fontsize=18 , loc='upper right' )
        #axarr[0].set_ylim(0.95,1.05) 
        axarr[1].set_xlim(-1 , 1)


['0.6', '0.8', '1', '1.2', '1.4']
---------------------------------------------------------------------------
IOError                                   Traceback (most recent call last)
<ipython-input-45-e575c6d08384> in <module>()
      9 ini =26
     10 for z in range(size(zvec2)):
---> 11         w , costi = readCost(ini, z)
     12         #w , costi = readCost(initials[ini], z)
     13         axarr[0].semilogx(w[w>0],costi[w>0]/costi[w<0][::-1])

<ipython-input-11-e4e90ccb988b> in readCost(ini, Ztrick)
    159         output = txt+repr(ini)+'/output_DMNRG_zEQ'+zvec2[Ztrick]+'_Dir'+repr(ini)+'.txt'
    160         #output = txt+repr(ini)+'/output_NRG_Code_DM_NRG_Dir'+repr(ini)+'.txt'
--> 161         infile = open(os.path.abspath(output), 'r')
    162         text = infile.readlines()
    163         w= []

IOError: [Errno 2] No such file or directory: '/home/cifucito/nrgcode/TwoChNRG/Data/run26/output_DMNRG_zEQ0.6_Dir26.txt'

In [35]:
fig, axMaj = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
fig2, axMaj2 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
vec0 = []
initials = array([70,71,72,73,74,75,76,77])
initials = array([57,58,59,60])
#initials = array([30,31,32,33,34])
#initials = array([35,36,37,38,39])
#initials = initials[::2]
for j in initials:
        #Read file
    vec , vec2 = readfile(j, 4)
    #axMaj.plot(vec2,vec, label= '$\epsilon_{d2} ='+str(-0.25+(j-initials[0])*0.05 )+'$') 
    axMaj.plot(vec2,vec, label= '$\epsilon_{d2} ='+str(0.005+(j-initials[0])*0.005 )+'$')
    vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
    vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]
    vec0.append(0.5*(vec0mas + vec0menos))


t2 = -0.25 + double(initials-initials[0]) * 0.05
t2 = 0.005+ double(initials-initials[0]) * 0.005
#t2[-1] = -0.06 ; t2[-2] = -0.08 
axMaj2.plot(t2,vec0, 'o' , label= title + '$(0)$')
  

axMaj.legend(fontsize=20 )
axMaj.set_xlim([-0.1,0.1])    
axMaj.set_xlabel('$\omega$',fontsize=25)
axMaj.set_ylabel('$\\rho_M$',fontsize=25)

axMaj2.set_xlabel('$\epsilon_{d2}$',fontsize=25)
axMaj2.set_ylabel('$\\rho_M(0)$',fontsize=25)
#axMaj2.legend(fontsize=20 , loc='upper left' )
axMaj2.set_xlim([0.005,0.02])
#axMaj2.set_ylim([0.0,4.5])


---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-35-b7a0cd8f92f3> in <module>()
     20 t2 = 0.005+ double(initials-initials[0]) * 0.005
     21 #t2[-1] = -0.06 ; t2[-2] = -0.08
---> 22 axMaj2.plot(t2,vec0, 'o' , label= title + '$(0)$')
     23 
     24 

TypeError: unsupported operand type(s) for +: 'function' and 'str'

In [ ]:
3/2

In [ ]:
L = array([5 , 2 , 4 ,10])

L[L >3]

In [ ]:
double(range(0,6)) * 0.05

In [ ]:
f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
f2, axarr2 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
f3, axarr3 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))



initials = array([40,42,43,44])
#initials = array([40,41,42,43 ,44,45])
initials = array([61,62,63,64,65])
#initials = initials[::2]
initials = array([56,57,58,59,60])
initials = array([139])


print initials
for i in range(4):
    vec0 = []
    #for j in range(40,46):
    for j in initials:
        #Read file
        vec , vec2 = readfile(j, i)
        #get 0 value of the DOS
        if j%5 == 1 : vec = vec[::-1];
        vec = array(vec) ; vec2 = array(vec2) 
        vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
        vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]   
        vec0.append(0.5*(vec0mas + vec0menos))
        axarr[i/2 , i%2 ].plot(vec2,vec, label= '$t_{2} ='+str(0.00+(j-initials[0])*0.005 )+'$') 
    
    
    axarr[i/2 , i%2 ].legend(fontsize=20)
    axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
    
    spin = '\uparrow' 
    if i%2 == 1 : spin = '\downarrow';  
    title = 'Dot $'+str(i/2 +1 )+'$ , Spin $'+ spin +'$'
    axarr[i/2 , i%2 ].set_title(title)    
    #axarr[i/2 , i%2 ].set_ylim([-0.0,9.0])
    if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\omega$',fontsize=25)
    if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('Spectral Density',fontsize=25)
        
    t2 = -0.0 + double(initials-initials[0]) * 0.005    
    axarr2.plot(t2,vec0 ,'o', label= title)  
    
    if i%2 == 0 : 
        new0 = vec0
    else: 
        print array(vec0)/array(new0) 
        axarr3.plot(t2,array(vec0)/array(new0) ,'o', label= title)
    
axarr2.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr2.set_ylabel('$\pi \Gamma \\rho(0)$',fontsize=25)
axarr2.legend(fontsize=20)
axarr2.set_xlim([-0.0,0.02])
axarr2.set_ylim([0.0,4.5])

axarr3.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr3.set_ylabel('$\\rho_\uparrow(0)/\\rho_\downarrow(0)$',fontsize=25)        
axarr3.legend(fontsize=20)        
axarr3.set_xlim([-0.0,0.02])
axarr3.set_ylim([0.0,1.8])

In [ ]:


In [8]:


In [8]:


In [10]:
f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
f2, axarr2 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))
f3, axarr3 = subplots(1, 1, sharex='col', sharey='row', figsize=(8, 5))



initials = array([40,42,43,44])
#initials = array([40,41,42,43 ,44,45])
initials = array([61,62,63,64,65])
#initials = initials[::2]
initials = array([35,36,37,38,39])
initials = initials[::2]
print initials
for i in range(4):
    vec0 = []
    #for j in range(40,46):
    for j in initials:
        #Read file
        vec , vec2 = readfile(j, i)
        
        #get 0 value of the DOS
        if j%20 == 1 : vec = vec[::-1];
        vec = array(vec) ; vec2 = array(vec2) 
        vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
        vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]   
        vec0.append(0.5*(vec0mas + vec0menos))
       
        axarr[i/2 , i%2 ].plot(vec2,vec, label= '$t_{2} ='+str(0.00+(j-initials[0])*0.005 )+'$') 
    
    
    axarr[i/2 , i%2 ].legend(fontsize=20)
    axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
    spin = '\uparrow' 
    if i%2 == 1 : spin = '\downarrow';  
    title = 'Dot $'+str(i/2 +1 )+'$ , Spin $'+ spin +'$'
    print title
    axarr[i/2 , i%2 ].set_title(title)    
    #axarr[i/2 , i%2 ].set_ylim([-0.0,9.0])
    if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\omega$',fontsize=25)
    if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('Spectral Density',fontsize=25)
        
    t2 = -0.0 + double(initials-initials[0]) * 0.005    
    axarr2.plot(t2,vec0 , label= title)  
    
    if i%2 == 0 : 
        new0 = vec0
    else: 
        print array(vec0)/array(new0) 
        axarr3.plot(t2,array(vec0)/array(new0) , label= title)
    
axarr2.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr2.set_ylabel('$\pi \Gamma \\rho(0)$',fontsize=25)
axarr2.legend(fontsize=20)
axarr2.set_xlim([-0.0,0.02])
axarr2.set_ylim([0.0,4.5])

axarr3.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr3.set_ylabel('$\\rho_\uparrow(0)/\\rho_\downarrow(0)$',fontsize=25)        
axarr3.legend(fontsize=20)        
axarr3.set_xlim([-0.0,0.02])
axarr3.set_ylim([0.0,1.8])


[35 37 39]
Dot $1$ , Spin $\uparrow$
Dot $1$ , Spin $\downarrow$
[ 0.50000003  0.50000001  0.5       ]
Dot $2$ , Spin $\uparrow$
Dot $2$ , Spin $\downarrow$
[ 0.50000026  0.50000004  0.5       ]
Out[10]:
(0.0, 1.8)

In [8]:


In [8]:


In [199]:
f, axarr = subplots(2, 2, sharex='col', sharey='row', figsize=(15, 10))
f2, axarr2 = subplots(1, 1, sharex='col', sharey='row', figsize=(15, 10))
f3, axarr3 = subplots(1, 1, sharex='col', sharey='row', figsize=(15, 10))

initials = array([30,31,32,33,34])
#initials = array([80,81,82,83,84])
initials = array([160,56])
#initials = initials[::2]
#initials = initials[4:5]
print initials
for i in range(4):
    vec0 = []
    #for j in range(40,46):
    for j in initials:
    #Read file
        vec , vec2 = readfile(j, i)
        #get 0 value of the DOS
        
        if j%20 == 2 : vec = vec[::-1];
        vec = array(vec) ; vec2 = array(vec2) 
        vec0mas = vec[ vec2 > 0] ; vec0mas = vec0mas[0]
        vec0menos = vec[ vec2 < 0] ; vec0menos = vec0menos[-1]
        vec0.append(0.5*(vec0mas + vec0menos))
        #print vec0mas , vec0menos , 0.5*(vec0mas + vec0menos)
        #maxima = vec[r_[True, vec[1:] > vec[:-1]] & r_[vec[:-1] > vec[1:], True]]
        #maxima = sort(maxima)
        
       
        axarr[i/2 , i%2 ].plot(vec2,vec,label= '$t_{2} ='+str(0.00+(j-initials[0])*0.005 )+'$') 
    print size(vec0)
    
    
    axarr[i/2 , i%2 ].legend(fontsize=20)
    axarr[i/2 , i%2 ].set_xlim([-0.5,0.5])
    spin = '\uparrow' 
    if i%2 == 1 : spin = '\downarrow';  
    title = 'Dot $'+str(i/2 +1 )+'$ , Spin $'+ spin +'$'
    print title
    axarr[i/2 , i%2 ].set_title(title)    
    #axarr[i/2 , i%2 ].set_ylim([-0.0,9.0])
    if i>1 : axarr[i/2 , i%2 ].set_xlabel('$\omega$',fontsize=25)
    if i%2 == 0 : axarr[i/2 , i%2 ].set_ylabel('Spectral Density',fontsize=25)
        
    t2 = -0.0 + double(initials-initials[0]) * 0.005    
    axarr2.plot(t2,vec0 , label= title)
    

        
    t2 = -0.0 + double(initials-initials[0]) * 0.005   
    
    if i%2 == 0 : 
        new0 = vec0
    else: 
        print array(vec0)/array(new0) 
        axarr3.plot(t2,array(new0)/array(vec0) , label= title)
    
axarr2.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr2.set_ylabel('$\pi \Gamma \\rho(0)$',fontsize=25)
axarr2.legend(fontsize=20)
axarr2.set_xlim([-0.0,0.02])
axarr2.set_ylim([0.0,4.5])

axarr3.set_xlabel('$t_1=t_{2}$',fontsize=25)
axarr3.set_ylabel('$\\rho_\uparrow(0)/\\rho_\downarrow(0)$',fontsize=25)        
axarr3.legend(fontsize=20)        
axarr3.set_xlim([-0.0,0.02])
#axarr3.set_ylim([0.0,1.0])


[160  56]
2
Dot $1$ , Spin $\uparrow$
2
Dot $1$ , Spin $\downarrow$
[ 1.  1.]
2
Dot $2$ , Spin $\uparrow$
2
Dot $2$ , Spin $\downarrow$
[ 1.  1.]
Out[199]:
(-0.0, 0.02)

In [ ]:
initials = array([109,119,124])
    for i in range(4):
        vec0 = []
        #for j in range(40,46):
        for j in initials:
            param = t2[where(initials == j)] ; param =  param[0]
            #Read file
            vec , vec2 = readfile(j, i)
            #get 0 value of the DOS