In [4]:

    

#from __future__ import division 
%pylab 
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


def readfile(ini):
        output = txt+repr(ini)+'/EntropyChain1Ch_NupPdn_zEQ1.dat'
        output = txt+repr(ini)+'/EntropyImp1ChNupPdn_Majorana_zEQ1.dat'
        infile = open(os.path.abspath(output), 'r')
        text = infile.readlines()
        
        #Data managemente
        Temp = [];         Tot = [];  Chain = [];  Imp = [];  
        for x in text:
            #print(list(x))
            a =x.split(' ') ; a[3]=a[3].split('\n')[0] 
            Temp.append(float(a[0]))
            Tot.append(float(a[1]))
            Chain.append(float(a[2]))
            Imp.append(float(a[3]))
        return array(Temp) , array(Tot) , array(Chain) , array(Imp)


    

    




Using matplotlib backend: TkAgg
Populating the interactive namespace from numpy and matplotlib






    




/usr/lib/python2.7/dist-packages/IPython/core/magics/pylab.py:161: UserWarning: pylab import has clobbered these variables: ['griddata']
`%matplotlib` prevents importing * from pylab and numpy
  "\n`%matplotlib` prevents importing * from pylab and numpy"

In [5]:

    

#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])   
    
    
    
#read files for folder ini and dot 
#return
def readfile2(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


    

In [6]:

    

ini = 152; line = 21
Temp , Tot , Chain , Imp = readfile(ini)
param = readfile_Param(ini, line)
#semilogx(Temp , Tot/log(2) , label = 'Tot'+ ', \eta_1 = '+ str(param) )
#semilogx(Temp , Chain/log(2) , label = 'Chain' + ', \eta_1 = '+ str(param) )
semilogx(Temp , Imp/log(2) , label = 'Imp'+  ',\eta_1 = '+ str(param) )
legend(fontsize=18,loc='lower right')
xlabel('Temperature')
ylabel('$\\times \log(2)$')


    

    
Out[6]:





Text(0,0.5,u'$\\times \\log(2)$')

In [13]:

    

legend(fontsize=23,loc='lower right')
#xlabel('Temperature')
#ylabel('$\\times \log(2)$')


    

    
Out[13]:





<matplotlib.legend.Legend at 0x7f6efe516b90>

In [13]:

    

ini = 150; line = 21
Temp , Tot , Chain , Imp = readfile(ini)
Rho , w  =readfile2(3, 0)
plot(range(size(w[w>0])),log(Rho[w>0]) )
plot(range(size(Imp[::-4])),Imp[::-4] )


    

    
Out[13]:





[<matplotlib.lines.Line2D at 0x7f6f770930d0>]

In [57]:

    

fill_between(linspace(min(params)/2 ,max(params)*2, 10), 0.002, 0.01,  facecolor='yellow')


    

    
Out[57]:





<matplotlib.lines.Line2D at 0x7f72cb7d8550>

In [7]:

    

def LogTransform(w , mult):
    logw = log10(w)
    return logw


    

In [8]:

    

from mpl_toolkits.axes_grid1 import make_axes_locatable


def LogMajSig( 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) 
    Params = linspace(min(t2) , max(t2),500) ;     SQ = 8*Params**2  
    Params= mult * LogTransform(Params , 1) ;     SQ = LogTransform(SQ , 1)
    
    KondoLims = array([0.002,0.01]) ; KondoLims = LogTransform(KondoLims , 1)
    ax.axvline(x=KondoLims[0], color = 'red', label = '$t_{1,2} =0$ Kondo Peaks')
    ax.axvline(x=KondoLims[1],color = 'red')
    #fill_between(linspace(KondoLims[1] ,KondoLims[2], 10),-4,-3, facecolor='green', label='Kondo Satellites')
    
    
    maxt2 = max(t2) ; mint2 = min(t2); t2 = []
    print maxt2 -mint2 
    
    for ini in initials:
        param =  readfile_Param(ini, line)
        param = mult * LogTransform(param , 1)
        t2.append(param)
        
        
        rho , w = readfile2(ini, dot) ;  rho1, w1 = readfile2(ini+100, dot)
        rho = 0.5*(rho +rho1) ; rho = rho[w>0] ; w = w[w>0]
        
        rhoD , w = readfile2(ini, dot+1) ;  rhoD1, w1 = readfile2(ini+100, dot+1)
        rhoD = 0.5*(rhoD +rhoD1) ; rhoD = rhoD[w>0] ; w = w[w>0]
        rho = rho/rhoD
        
        
        maxi = max(-log10(w[w>0]))
        logw = LogTransform(w, 1) ; 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.5, vmax=maxplot,cmap="gist_stern", 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]  
    ax.plot(SQ , Params, 'b-',label = '$4t^2/U$')
    ax.legend(fontsize=20, loc='upper left')
    
    if PrintAxis :
        param_ticks = array([10**(-4), 10**(-3), 10**(-2)])
        Y_logticks = mult * LogTransform(param_ticks , 1)
        yticks(Y_logticks , ('$10^{-4}$','$10^{-3}$', '$10^{-2}$')) 
        #print 'ParamTicks =' , paramTicks_scaled , paramTicks
        
        ticks = array([10**(-8) , 10**(-6) ,10**(-4) , 10**(-2),1])
        logticks = LogTransform(ticks , 1)
        xticks(logticks , ('$10^{-8}$','$10^{-6}$','$10^{-4}$','$10^{-2}$', '$1$'))     
        cbar = f2.colorbar(im)
        cbar.set_ticks([0.5, 1,1.5, 2])
        cbar.set_ticklabels(['$0.5$','$1$' ,'$1.5$', '$2$'])
    #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 logticks  , Y_logticks


    

In [ ]:

    




    

In [9]:

    

def LogTransform2(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 [44]:

    

def LogMajt1t2( 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) 
    Params = linspace(min(log10(t2)) , max(log10(t2)),20) ; Params2 = 10**Params ;     SQ = 8*Params2**2  
     #LogTransform(SQ , 1)
    
    KondoLims = array([0.002,0.01]) ; KondoLims = LogTransform(KondoLims , 1)
#     ax.axvline(x=KondoLims[0], color = 'red', label = '$t_{1,2} =0$ Kondo Peaks')
#     ax.axvline(x=KondoLims[1],color = 'red')
    #fill_between(linspace(KondoLims[1] ,KondoLims[2], 10),-4,-3, facecolor='green', label='Kondo Satellites')
    
    
    maxt2 = max(t2) ; mint2 = min(t2); t2 = []
    print maxt2 -mint2 
    
    for ini in initials:
        param =  readfile_Param(ini, line)
        param = mult * LogTransform(param , 1)
        t2.append(param)
        
        
        rho , w = readfile2(ini, dot) ;  rho1, w1 = readfile2(ini+100, dot)
        rho = 0.5*(rho +rho1) ;# rho = rho[w>0] ; w = w[w>0]
#         rho = rho*pi*Gamma ; w = w/Gamma
        
        
#         rhoD , w = readfile2(ini, dot+1) ;  rhoD1, w1 = readfile2(ini+100, dot+1)
#         rhoD = 0.5*(rhoD +rhoD1) ; rhoD = rhoD[w>0] ; w = w[w>0]
 
        maxi = max(-log10(w[w>0]))
        logw = LogTransform2(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.5, vmax=maxplot,cmap="gist_stern", 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]  
#     ax.plot(SQ , Params, 'b-',label = '$4t^2/U$')
#         ax.legend(fontsize=20, loc='upper left')
    Params= mult * Params ;     SQ = LogTransform2(SQ, maxi)
    ax.plot(SQ , Params, 'o' , color = 'g' ,label = '$4t^2/U$')
    ax.plot(-SQ , Params, 'o' , color = 'g' ,label = '$4t^2/U$')
    ax.legend(fontsize=20, loc='lower left')
    
    
    if PrintAxis :
        param_ticks = array([10**(-4), 10**(-3), 10**(-2)])
        Y_logticks = mult * LogTransform(param_ticks , 1)
        yticks(Y_logticks , ('$10^{-4}$','$10^{-3}$', '$10^{-2}$')) 
    
        ticks = array([-1,-10**(-3) ,-10**(-6) ,10**(-maxi), 10**(-6) ,  10**(-3),1])
        logticks = LogTransform2(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


    

In [45]:

    

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

#initials = array([180,181,182,183,184,190 ,191,192,193,194,195,186,187,188,189])
initials = range(811,820)
f2, ax = subplots(1, 1, sharex='col', sharey='row', figsize=(7, 6)) 
spin = '\downarrow' if dot%2 == 1  else '\uparrow';
ax.set_title('$\\rho_'+ spin +'$')
ax.set_xlabel('$\omega / \Gamma$')
ax.set_ylabel('$t_1/D =t_2 / \Gamma$')

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 
logticks  , Y_logticks  = LogMajt1t2( initials , dot , line , f2, ax, Method, 2.5, True , 2.5)


    

    




0.09995

In [ ]:

    




    

In [ ]:

    




    

In [29]:

    

#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( Plus ,   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)/Gamma ;  t2.append(param)
    print t2    
    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 ]
        
        maxi , im = LogMajt1t2( 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 = LogTransform2(ticks , maxi)
        wTicksLabels = ['$-1$','$-10^{-3}$', '$-10^{-6}$','$0$', '$10^{-6}$','$10^{-3}$', '$1$']
            
        param_ticks = array([10**(-4), 10**(-3), 10**(-2)])
        Y_logticks = mult * LogTransform(param_ticks , 1)
        yticks(Y_logticks , ('$10^{-4}$','$10^{-3}$', '$10^{-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=['$10^{-4}$','$10^{-3}$', '$10^{-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/\Gamma'+ '$')
        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)

    

Method = 'linear' #; Method = 'nearest';  
line = 13
xtitle = [None]*23
xtitle[4] = '\epsilon_{d1}/\Gamma' ; xtitle[12] = 't_1/\Gamma = t_2/\Gamma' ; xtitle[12] = 't_1/\Gamma' ;
xtitle[16] = '\epsilon_{d2}/\Gamma'; xtitle[15] = '\epsilon_{d2}/\Gamma = \\frac{-U_2}{\Gamma}';xtitle[13] = 't_2/\Gamma'
xtitle[17] = 't_{dots}/\Gamma ' ; xtitle[18] = '\epsilon_m/\Gamma ' ;xtitle[19] = '\phi_1' ;  xtitle[20] = '\phi_2' ;xtitle[21] = '\eta_1' ;
xtitle[22] = '\eta_1' 
rcParams.update({'font.size': 30})
initials = array([450,451,452,453,454,455,456,457,458,459,816,817,818])
initials = range(811,819)



color4Dots(True ,1,initials ,  line , wlims, xtitle[line], Method ,0.25,1.5)


    

    




[0.0017687156648071572, 0.0035374313296143143, 0.01768715664807157, 0.03537431329614314, 0.1768715664807157, 0.3537431329614314, 1.0, 1.7687156648071571]
1.76694694914 [ 0.0015015   0.01501502  0.15015015  0.84892192]
0.04995
[ 0.0015015   0.42521171  0.84892192] [ 0.00176872  0.50088436  1.        ]
0.04995
[ 0.0015015   0.42521171  0.84892192] [ 0.00176872  0.50088436  1.        ]






    




Exception in Tkinter callback
Traceback (most recent call last):
  File "/usr/lib/python2.7/lib-tk/Tkinter.py", line 1544, in __call__
    return self.func(*args)
  File "/usr/lib/python2.7/lib-tk/Tkinter.py", line 595, in callit
    func(*args)
  File "/usr/lib/python2.7/dist-packages/matplotlib/backends/backend_tkagg.py", line 323, in idle_draw
    self.draw()
  File "/usr/lib/python2.7/dist-packages/matplotlib/backends/backend_tkagg.py", line 304, in draw
    FigureCanvasAgg.draw(self)
  File "/usr/lib/python2.7/dist-packages/matplotlib/backends/backend_agg.py", line 430, in draw
    self.figure.draw(self.renderer)
  File "/usr/lib/python2.7/dist-packages/matplotlib/artist.py", line 55, in draw_wrapper
    return draw(artist, renderer, *args, **kwargs)
  File "/usr/lib/python2.7/dist-packages/matplotlib/figure.py", line 1299, in draw
    renderer, self, artists, self.suppressComposite)
  File "/usr/lib/python2.7/dist-packages/matplotlib/image.py", line 138, in _draw_list_compositing_images
    a.draw(renderer)
  File "/usr/lib/python2.7/dist-packages/matplotlib/artist.py", line 55, in draw_wrapper
    return draw(artist, renderer, *args, **kwargs)
  File "/usr/lib/python2.7/dist-packages/matplotlib/axes/_base.py", line 2437, in draw
    mimage._draw_list_compositing_images(renderer, self, artists)
  File "/usr/lib/python2.7/dist-packages/matplotlib/image.py", line 138, in _draw_list_compositing_images
    a.draw(renderer)
  File "/usr/lib/python2.7/dist-packages/matplotlib/artist.py", line 55, in draw_wrapper
    return draw(artist, renderer, *args, **kwargs)
  File "/usr/lib/python2.7/dist-packages/matplotlib/image.py", line 555, in draw
    renderer, renderer.get_image_magnification())
  File "/usr/lib/python2.7/dist-packages/matplotlib/image.py", line 782, in make_image
    unsampled=unsampled)
  File "/usr/lib/python2.7/dist-packages/matplotlib/image.py", line 451, in _make_image
    output = self.norm(np.ma.masked_array(A_resampled, out_mask))
  File "/usr/lib/python2.7/dist-packages/matplotlib/colors.py", line 928, in __call__
    raise ValueError("minvalue must be less than or equal to maxvalue")
ValueError: minvalue must be less than or equal to maxvalue

In [ ]:

    




    

In [ ]:

    




    

In [ ]:

    




    

In [7]:

    

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

#initials = array([180,181,182,183,184,190 ,191,192,193,194,195,186,187,188,189])
initials = range(811,819)
f2, ax = subplots(1, 1, sharex='col', sharey='row', figsize=(10, 10))
ax.set_title('Majorana Signature $\\rho_\uparrow / \\rho_\downarrow$', fontsize=22)
ax.set_xlabel('$\omega / D$')
ax.set_ylabel('$t_1/D =t_2 / D$')

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 
logticks  , Y_logticks  = LogMajSig( initials , dot , line , f2, ax, Method, 2, True , 2.5)


#print maxrho  , maxplot


    

    




0.04995

In [8]:

    

ax.legend(fontsize=19, loc='upper left')
ax.set_title('Majorana Signature $\\rho_\uparrow / \\rho_\downarrow$', fontsize=22)
xticks(logticks , ('$10^{-8}$','$10^{-6}$','$10^{-4}$','$10^{-2}$', '$1$'))
yticks(Y_logticks , ('$10^{-4}$','$10^{-3}$', '$10^{-2}$'))


    

    
Out[8]:





([<matplotlib.axis.YTick at 0x7f94b11ebe90>,
  <matplotlib.axis.YTick at 0x7f94b11e2e90>,
  <matplotlib.axis.YTick at 0x7f94b1e27e50>],
 <a list of 3 Text yticklabel objects>)

In [74]:

    

def GetTemp(ini , dot , trigger):
    dot = (dot-1)*2 
    rho , w = readfile2(ini,dot) ; rho1 , w = readfile2(ini+100,dot) ; rho = 0.5*(rho +rho1)
    rhoD , w = readfile2(ini,dot+1) ; rhoD1 , w = readfile2(ini+100,dot+1) ; rhoD = 0.5*(rhoD +rhoD1)
    Temp = w[rho/rhoD>trigger] 
    #semilogx(w,rho,w , rhoD)
    return Temp[-1]
    #semilogx(w,rho/rhoD)
    
ini = range(811,819)
Temp1 = [] ;  Temp2 = [] ; params = []
for i in ini:
    Temp1.append(GetTemp(i, 1 , 1.3) )
    Temp2.append(GetTemp(i, 1 , 1.95) )#/Gamma)
    params.append(readfile_Param(i,12))
    
params = array(params)    



x = linspace(min(params)/2 ,max(params)*2, 200)

#loglog(tnew,f(tnew),'-' )
loglog(params, Temp1 ,  'bo' , label='Limits of Majorana Satellites' )
loglog(params, Temp2 ,  'bo' )#, label='Minimum Majorana Satellites' )
loglog(x , 8*(x**2), 'r-' , label = '$4t^2/U$')
#axhline(y=0.002, color='r', linestyle='-' , )
#axhline(y=0.01, color='r', linestyle='-')
fill_between(linspace(min(params)/2 ,max(params)*2, 10), 0.002, 0.01,  facecolor='green', label='Kondo Satellites')
fill_between(params, Temp1, Temp2,  facecolor='silver' , label=' Majorana Satellite \\rho_\uparrow / \\rho_\downarrow >1.5')
fill_between(params, min(Temp2), Temp2,  facecolor='gold' , label='Majorana Signature \\rho_\uparrow / \\rho_\downarrow =2')
#fill_between(params, 0., Temp1,  facecolor='red' , label='Majorana Peaks')
xlim(min(params) , max(params)) 
ylim(min(Temp2) , max(Temp1)) 
#for p in range(len(params)):
    ##axvline(x=params[p], ymin= 2, ymax = 30)
    #ymax=Temp1[p])
#    print Temp2[p] , Temp1[p]

#rho , w = readfile2(450,dot) ; rho1 , w = readfile2(450+100,dot) ; rho = 0.5*(rho +rho1)
legend(fontsize=18, loc='lower right')
xlabel('$t_1/D = t_2/D$')
ylabel('T_M / D)')


    

    
Out[74]:





<matplotlib.text.Text at 0x7fed813872d0>

In [75]:

    

legend(fontsize=18, loc='lower right')


    

    
Out[75]:





<matplotlib.legend.Legend at 0x7fed80cd61d0>

In [172]:

    

from scipy import interpolate
from scipy.optimize import curve_fit
def func(x,a,b):
     return a*x +b 
    
tnew =  linspace(min(params) ,max(params), 100) #; f = interpolate.interp1d(params, Temp1 , kind='quadratic')
popt, pcov = curve_fit(func, log(params), log(Temp1))   
data=popt[1]*tnew**popt[0]
#loglog(tnew ,data)
semilogx(tnew,data)
print tnew , data


    

    




[  5.00000000e-05   5.54545455e-04   1.05909091e-03   1.56363636e-03
   2.06818182e-03   2.57272727e-03   3.07727273e-03   3.58181818e-03
   4.08636364e-03   4.59090909e-03   5.09545455e-03   5.60000000e-03
   6.10454545e-03   6.60909091e-03   7.11363636e-03   7.61818182e-03
   8.12272727e-03   8.62727273e-03   9.13181818e-03   9.63636364e-03
   1.01409091e-02   1.06454545e-02   1.11500000e-02   1.16545455e-02
   1.21590909e-02   1.26636364e-02   1.31681818e-02   1.36727273e-02
   1.41772727e-02   1.46818182e-02   1.51863636e-02   1.56909091e-02
   1.61954545e-02   1.67000000e-02   1.72045455e-02   1.77090909e-02
   1.82136364e-02   1.87181818e-02   1.92227273e-02   1.97272727e-02
   2.02318182e-02   2.07363636e-02   2.12409091e-02   2.17454545e-02
   2.22500000e-02   2.27545455e-02   2.32590909e-02   2.37636364e-02
   2.42681818e-02   2.47727273e-02   2.52772727e-02   2.57818182e-02
   2.62863636e-02   2.67909091e-02   2.72954545e-02   2.78000000e-02
   2.83045455e-02   2.88090909e-02   2.93136364e-02   2.98181818e-02
   3.03227273e-02   3.08272727e-02   3.13318182e-02   3.18363636e-02
   3.23409091e-02   3.28454545e-02   3.33500000e-02   3.38545455e-02
   3.43590909e-02   3.48636364e-02   3.53681818e-02   3.58727273e-02
   3.63772727e-02   3.68818182e-02   3.73863636e-02   3.78909091e-02
   3.83954545e-02   3.89000000e-02   3.94045455e-02   3.99090909e-02
   4.04136364e-02   4.09181818e-02   4.14227273e-02   4.19272727e-02
   4.24318182e-02   4.29363636e-02   4.34409091e-02   4.39454545e-02
   4.44500000e-02   4.49545455e-02   4.54590909e-02   4.59636364e-02
   4.64681818e-02   4.69727273e-02   4.74772727e-02   4.79818182e-02
   4.84863636e-02   4.89909091e-02   4.94954545e-02   5.00000000e-02] [  4.67185834e-08   3.18754177e-06   9.92240995e-06   1.96596375e-05
   3.21167153e-05   4.71105735e-05   6.45078112e-05   8.42045944e-05
   1.06116572e-04   1.30173125e-04   1.56313788e-04   1.84485886e-04
   2.14642891e-04   2.46743235e-04   2.80749431e-04   3.16627405e-04
   3.54345971e-04   3.93876417e-04   4.35192175e-04   4.78268543e-04
   5.23082467e-04   5.69612351e-04   6.17837894e-04   6.67739965e-04
   7.19300479e-04   7.72502302e-04   8.27329163e-04   8.83765580e-04
   9.41796792e-04   1.00140870e-03   1.06258783e-03   1.12532124e-03
   1.18959655e-03   1.25540185e-03   1.32272567e-03   1.39155699e-03
   1.46188517e-03   1.53369994e-03   1.60699138e-03   1.68174990e-03
   1.75796622e-03   1.83563135e-03   1.91473657e-03   1.99527343e-03
   2.07723371e-03   2.16060944e-03   2.24539286e-03   2.33157644e-03
   2.41915282e-03   2.50811487e-03   2.59845561e-03   2.69016825e-03
   2.78324616e-03   2.87768288e-03   2.97347209e-03   3.07060763e-03
   3.16908347e-03   3.26889371e-03   3.37003260e-03   3.47249449e-03
   3.57627386e-03   3.68136532e-03   3.78776356e-03   3.89546340e-03
   4.00445976e-03   4.11474765e-03   4.22632219e-03   4.33917857e-03
   4.45331210e-03   4.56871815e-03   4.68539218e-03   4.80332975e-03
   4.92252647e-03   5.04297804e-03   5.16468024e-03   5.28762891e-03
   5.41181996e-03   5.53724938e-03   5.66391321e-03   5.79180757e-03
   5.92092861e-03   6.05127257e-03   6.18283575e-03   6.31561448e-03
   6.44960517e-03   6.58480427e-03   6.72120829e-03   6.85881378e-03
   6.99761735e-03   7.13761566e-03   7.27880540e-03   7.42118332e-03
   7.56474621e-03   7.70949092e-03   7.85541431e-03   8.00251331e-03
   8.15078488e-03   8.30022601e-03   8.45083374e-03   8.60260516e-03]

In [74]:

    

GetTemp(810,1)


    

In [72]:

    

dot = 0
i = 115
rho , w = readfile2(i,dot) ; rho1 , w = readfile2(i+100,dot) ; rho = 0.5*(rho +rho1)
rhoD , w = readfile2(i,dot+1) ; rhoD1 , w = readfile2(i+100,dot+1) ; rhoD = 0.5*(rhoD +rhoD1)
peaks = r_[True, rho[1:] < rho[:-1]] & r_[rho[:-1] < rho[1:], True]
KondoTemp = w[rho>(0.5*max(rho))] 
print KondoTemp[-1]
semilogx(w, rho, w, rhoD)
#w[peaks]


    

    




0.000175283469751






    
Out[72]:





[<matplotlib.lines.Line2D at 0x7f72c9d39510>,
 <matplotlib.lines.Line2D at 0x7f72c9956d50>]

In [78]:

    

dot = 0
#readfile2(450,dot)


    

In [73]:

    

ini = 450  ; dot = 0 
    rho , w = readfile2(ini,dot) ; rho1 , w = readfile2(ini+100,dot) ; rho = 0.5*(rho +rho1)
    rhoD , w = readfile2(ini,dot+1) ; rhoD1 , w = readfile2(ini+100,dot+1) ; rhoD = 0.5*(rhoD +rhoD1)
    Temp = w[rho/rhoD>1.5] 
    semilogx(w, rho/rhoD , w , rho , w, rhoD)


    

    
Out[73]:





[<matplotlib.lines.Line2D at 0x7fed88d22e90>,
 <matplotlib.lines.Line2D at 0x7fed818160d0>,
 <matplotlib.lines.Line2D at 0x7fed81816790>]

In [ ]: