MinimalConductivity.m

Go to the documentation of this file.

%%    Minimal conductivity of graphene
%    Copyright (C) 2015 Peter Rakyta, Ph.D.
%
%    This program is free software: you can redistribute it and/or modify
%    it under the terms of the GNU General Public License as published by
%    the Free Software Foundation, either version 3 of the License, or
%    (at your option) any later version.
%
%    This program is distributed in the hope that it will be useful,
%    but WITHOUT ANY WARRANTY; without even the implied warranty of
%    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%    GNU General Public License for more details.
%
%    You should have received a copy of the GNU General Public License
%    along with this program.  If not, see http://www.gnu.org/licenses/.
%
%> @addtogroup Examples
%> @{
%> @file MinimalConductivity.m
%> @brief Calculates the minimal conductivity of the graphene sheet as a function of the aspect ratio.
%> @param filenum The identification number of the filenema for the exported data (default is 1).
%> @Available
%> EQuUs v4.4 or later
%> @expectedresult
%> @image html MinimalConductivity.jpg
%> @image latex MinimalConductivity.jpg
%> The red line represents the theoretical \f$ 2/\pi\;\sigma_0 \f$ limit of the minimal conductivity with \f$ \sigma_0 = e^2/\hbar \f$.
%> @}
%
%> @brief Calculates the minimal conductivity of the graphene sheet as a function of the aspect ratio.
%> @param filenum The identification number of the filenema for the exported data (default is 1).
function MinimalConductivity( filenum )

if ~exist('filenum', 'var')
    filenum = 1;
end

filename = mfilename('fullpath');
[directory, fncname] = fileparts( filename );

    % filename containing the input XML    
    inputXML = 'Basic_Input_zigzag_leads.xml';
    % Parsing the input file and creating data structures
    [Opt, param] = parseInput( fullfile( directory, inputXML ) );

    % filename containing the output XML
    outfilename = [fncname, '_',num2str( filenum )];
    % the output folder
    outputdir   = [];    

    % length of the junction (number of unit cells)
    height = 170;
    % The Fermi energy
    EF = 0.001; %In EV

    % class representing the two terminal ribbon structure
    cRibbon    = [];


    % creating 1D energy array
    % The minimal value of the widths
    Widthmin = 10;
    % spacing of the width values
    Widthdelta = 20;
    % number of width values
    Widthnum = 30;%100;%200;
    % The 1D energy array
    Widthvec = Widthmin:Widthdelta:Widthmin+Widthnum*Widthdelta;
    
    % preallocating the array for the calculated conductance values
    Conductivity = NaN(Widthnum+1,1); 
    % preallocating the array for the unitarity errors
    DeltaC      = NaN(Widthnum+1,1);
    % preallocating the array for the spect ratio (Width/lenth)
    aspectRatio = NaN(Widthnum+1,1); 

    % creating output directories
    setOutputDir();

    % Calculate the transport through the magnetic barrier
    CalculateTransport();
    
    % plot the calculated results
    PlotFunction();




%% CalculateTransport
%> @brief Calculates the conductance values
    function CalculateTransport()
        
        
        
        for idx = 1:length(Widthvec)  
            
            % the current width parameter of the junction
            width = Widthvec(idx);
            
            % creating the Ribbon class representing the twoterminal setup
            cRibbon = Ribbon('width', width, 'height', height, 'Opt', Opt, 'param', param, 'EF', EF );     
        
            disp(' ----------------------------------')
            disp('Plotting conductivity:')
            PlotFunction();        
        
        
            disp(' ')
            disp(' ')
            disp(' ----------------------------------')
            disp('Calculating Conductance for the given aspect ratio')           
        

            % Calculating the conductivity
            tic
            [Conductivity_tmp,aspectRatio_tmp,~,~,DeltaC_tmp] = cRibbon.Transport( 0, 'constant_channels', true );
            toc
            
            % storing the calculated data
            Conductivity(idx) = Conductivity_tmp;
            DeltaC(idx) = DeltaC_tmp;
            aspectRatio(idx) = aspectRatio_tmp;
            
            
            
            % exporting the calculated data
            save( fullfile(outputdir,[outfilename, '.mat']) );
            disp([ num2str(idx/length(Widthvec)*100),' % calculated of the conductance.'])
            disp( ' ' )
            
            
            
        end
                
    end



%% setOutputDir
%> @brief Sets output directory.
    function setOutputDir()
        resultsdir = [pwd, filesep, 'results'];
        mkdir(resultsdir );
        outputdir = resultsdir;        
    end



 %% PlotFunction
%> @brief Creates the plot
    function PlotFunction( )
        
        % creating figure in units of pixels
        figure1 = figure( 'Units', 'Pixels', 'Visible', 'off');
        
        % font size on the figure will be 16 points
        fontsize = 16;
          
        
        %Theoretical minimal conductivity
        theoretical_limit = 2/pi; 
        
        % determine points to be plotted
        indexek = logical( (abs(DeltaC) <= 0.5) & (aspectRatio~=0) & ~isnan(DeltaC) );
        x_lim = [0 max(aspectRatio(indexek))];
        
        if norm(x_lim) == 0
            x_lim = [0 1];
        end
        
        axes_cond = axes('Parent',figure1, ...
                'Visible', 'on',...
                'FontSize', fontsize,... 
                'xlim', x_lim,...   
                'Box', 'on',...
                'Units', 'Pixels', ...
                'FontName','Times New Roman');
        hold on; 
        
        % plot the data
        numerics = plot(aspectRatio(indexek), Conductivity(indexek), 'Linewidth', 2, 'color', [0 0 0], 'Parent', axes_cond);   
        try
            plot(aspectRatio(indexek), DeltaC(indexek), 'Linewidth', 1, 'color', [0 1 0], 'Parent', axes_cond);
        catch
            disp('No DeltaC data present')
        end
        plot( [0 max(aspectRatio(indexek))], ones(1,2)*theoretical_limit, 'r',  'Parent', axes_cond);
        
               
         
        
         % Create xlabel
         xlabel('W/L','FontSize', fontsize,'FontName','Times New Roman', 'Parent', axes_cond);

         % Create ylabel
         ylabel('\sigma/\sigma_0','FontSize', fontsize,'FontName','Times New Roman', 'Parent', axes_cond);
        
        
         % set the legends
         fig_legend = legend(axes_cond, {'conductivity', 'unitarity error'});%, 'FontSize', fontsize, 'FontName','Times New Roman')
         set(fig_legend, 'FontSize', fontsize, 'FontName','Times New Roman', 'Box', 'off', 'Location', 'NorthEast')
         
         % setting the position and margins of the plot, removing white
         % spaces for release dates greater than 2015
         ver = version('-release');
         if str2num(ver(1:4)) >= 2016
            figure1.PaperPositionMode = 'auto';
            fig_pos = figure1.PaperPosition;
            figure1.PaperSize = [fig_pos(3) fig_pos(4)]; 
        
            set(axes_cond, 'Position', get(axes_cond, 'OuterPosition') - get(axes_cond, 'TightInset') * [-1 0 1 0; 0 -1 0 1; 0 0 1 0; 0 0 0 1]);        
            Position_figure = get(axes_cond, 'OuterPosition');
            set(figure1, 'Position', Position_figure);
         end

        
         % export the figures
         print('-depsc2', fullfile(outputdir,[outfilename, '.eps']))
         print('-dpdf', fullfile(outputdir,[outfilename, '.pdf']))
         close(figure1)
        
        
    end




end