Graphene_QAD_eig.m

Go to the documentation of this file.

%    Calculates the eigenstates of quantum antidots - based on EQuUs v4.4
%    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 Graphene_QAD_eig.m
%> @brief Example to calculate the eigenstates of quantum antidots.
%> @param filenum The identification number of the filenema for the exported data (default is 1).
%> @Available
%> EQuUs v4.7 or later
%> @expectedresult
%> @image html Graphene_QAD_eig.jpg
%> @image latex Graphene_QAD_eig.jpg
%> @}
%
%> @brief Example to calculate the eigenstates of quantum antidots.
%> @param filenum The identification number of the filenema for the exported data (default is 1).
function Graphene_QAD_eig( 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   = [];      

     %> The width of the graphene strip
    width = 320;
     %> The length of the graphene strip
    height = 350;
     %> The center and radius of the antidot
    hole = struct( 'center', [width/2 height/2], 'radius', 90);  
    %> Magnetic field in Tesla    
    B = 40;
     % the output directory
     outputdir   = [];
    %> The calculated wave functions of the antidot
     hole_wavefnc = [];

     %> setting the output directory
     setOutputDir()
    
     %> An instance of class antidot
     cAntidot = antidot('width', width, 'height', height, 'Opt', Opt, 'param', param, 'B', B, 'hole', hole, ...
                        'rCC', 1.42, 'filenameOut', fullfile( outputdir, [outfilename, '.xml']) );
    
     %> calculate the wave functions
     hole_wavefnc = calc_TBwavefnc();


     % plot the wave functions if screen is on.
    plotWavefnc();


%% calc_TBwavefnc
% @brief Calculates TB wave functions
    function hole_wavefnc = calc_TBwavefnc()
        
          % set filterHole=1 to sort out the wave functions localized to the boundaries of the graphene strip       
        filterHole = 1;

          % expected energies (in eV) of the antidot states for different edges
        Energy_infmass = 0.1233529;
        Energy_smooth = 0.1102997;
        Energy_rough  = 0.1080343; 

          % calculate the eigenstate for infinite mass boundary condition
        hole_wavefnc.infmass = cAntidot.CalcWavefnc( Energy_infmass, 1, 'toPlot', 0, 'filterHole', filterHole, 'db', 1, 'infinitemass', 1, 'dotCalc', 0 );

          % calculate the eigenstate for smooth edge boundary condition
        hole_wavefnc.smoothedge = cAntidot.CalcWavefnc( Energy_smooth, 1, 'toPlot', 0, 'filterHole', filterHole, 'db', 1, 'infinitemass', 0, 'dotCalc', 0, 'smoothedge', 1 ); 

          % calculate the eigenstate for rough edge boundary condition
        hole_wavefnc.roughedge = cAntidot.CalcWavefnc( Energy_rough, 1, 'toPlot', 0, 'filterHole', filterHole, 'db', 1, 'infinitemass', 0, 'dotCalc', 0 );

          % export the calculetd wave functions
        save( [outputdir, filesep, outfilename,'.mat'], 'hole_wavefnc' );        
        
    end


%% plotWavefnc
% @brief plots the calculated wave functions
    function plotWavefnc()
    
        disp( 'Plotting wave functions' )   
        
        % if there is no display, the plot is skipped
        if usejava('jvm') && ~feature('ShowFigureWindows')
            disp('No display detected, exiting from plot.')
            return
        end
    
        % define custom colormap levels
        Colormap_levels = [1 1 0.989610493183136;0.910812675952911 0.910812675952911 0.86870539188385;0.868301570415497 0.868301570415497 0.812190175056458;0.825790464878082 0.825790464878082 0.755674958229065;0.809619069099426 0.809619069099426 0.734643757343292;0.793447613716125 0.793447613716125 0.713612616062164;0.777276217937469 0.777276217937469 0.692581415176392;0.761104762554169 0.761104762554169 0.671550273895264;0.744933366775513 0.744933366775513 0.650519073009491;0.725760400295258 0.725760400295258 0.626198768615723;0.706587433815002 0.706587433815002 0.601878464221954;0.687414407730103 0.687414407730103 0.577558159828186;0.668241441249847 0.668241441249847 0.553237855434418;0.650109708309174 0.650109708309174 0.530886232852936;0.6319779753685 0.6319779753685 0.508534550666809;0.613846242427826 0.613846242427826 0.486182928085327;0.595714509487152 0.595714509487152 0.463831305503845;0.57862401008606 0.57862401008606 0.443448305130005;0.561533570289612 0.561533570289612 0.423065304756165;0.544443130493164 0.544443130493164 0.402682334184647;0.527352631092072 0.527352631092072 0.382299333810806;0.505953729152679 0.505953729152679 0.35774689912796;0.484554797410965 0.484554797410965 0.333194434642792;0.463155895471573 0.463155895471573 0.308641999959946;0.450895965099335 0.446475297212601 0.298696875572205;0.438636004924774 0.429794728755951 0.288751751184464;0.426376074552536 0.413114130496979 0.278806626796722;0.414116144180298 0.396433532238007 0.268861502408981;0.402542054653168 0.381049633026123 0.259602248668671;0.390967965126038 0.36566573381424 0.250342965126038;0.379393875598907 0.350281864404678 0.241083696484566;0.367819786071777 0.334897965192795 0.231824427843094;0.3569315969944 0.320810824632645 0.223251044750214;0.346043407917023 0.306723684072495 0.214677661657333;0.335155189037323 0.292636543512344 0.206104278564453;0.324266999959946 0.278549402952194 0.197530895471573;0.314064651727676 0.265758991241455 0.189643383026123;0.303862333297729 0.252968549728394 0.181755855679512;0.29365998506546 0.240178138017654 0.173868328332901;0.283457636833191 0.227387711405754 0.165980815887451;0.273941189050674 0.215894013643265 0.158779174089432;0.264424741268158 0.204400300979614 0.151577517390251;0.254908263683319 0.192906603217125 0.144375860691071;0.245391815900803 0.181412905454636 0.137174218893051;0.236561223864555 0.171215921640396 0.130658447742462;0.227730631828308 0.161018952727318 0.124142669141293;0.218900039792061 0.15082198381424 0.117626890540123;0.210069447755814 0.140625 0.111111119389534;0.201924726366997 0.131724759936333 0.105281211435795;0.19378000497818 0.122824504971504 0.0994513034820557;0.185635298490524 0.113924264907837 0.0936214029788971;0.177490577101707 0.105024017393589 0.0877914950251579;0.170031726360321 0.0974204912781715 0.0826474577188492;0.162572875618935 0.0898169651627541 0.077503427863121;0.155114024877548 0.0822134390473366 0.0723593980073929;0.147655174136162 0.0746099129319191 0.0672153607010841;0.140882194042206 0.0683031156659126 0.0627572014927864;0.13410921394825 0.0619963146746159 0.0582990385591984;0.127336233854294 0.0556895136833191 0.0538408756256104;0.120563261210918 0.0493827164173126 0.0493827164173126;0.0904224440455437 0.0370370373129845 0.0370370373129845;0.0602816306054592 0.0246913582086563 0.0246913582086563;0.0301408153027296 0.0123456791043282 0.0123456791043282;0 0 0];
 
        % creating figure in units of pixels and utilizing the previous colormap
        figure1 = figure( 'Units', 'Pixels','Renderer','painters', 'Visible', 'off', ...
                    'Colormap', Colormap_levels);
        
        % font size on the figure in points
        fontsize = 14;
        
        % 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
            fig = gcf;
            fig.PaperPositionMode = 'auto';
            fig_pos = fig.PaperPosition;
            fig.PaperSize = [fig_pos(3) fig_pos(4)]; 
            fig_pos(4) = fig_pos(4)/3;        
        end

        
        y_lim = [200 650];
        delta = 1;
        
        
        %************** wave function infinite mass **********************
        
        axes_geo1 = axes( 'Parent',figure1,...
                'CLim',[0 1e-4],...
                'Visible', 'off',...
                'FontSize', fontsize,...    
                'ylim', y_lim,...   
                'Units', 'Pixels', ...
                'FontName','Times New Roman');
        hold on; 
        axis equal;
        
        rCC = cAntidot.rCC*1e-10; %Angstrom
        phi0 = cAntidot.h/cAntidot.qe; % magnetic flux quantum
        eta_B = 2*pi/phi0*(rCC)^2*B;                          
        homega = 2.97*sqrt(9/2*eta_B);
        
         % selecting data to plot
        Wavefnc = hole_wavefnc.infmass.Wavefnc{1};
        xcoords = hole_wavefnc.infmass.xcoords;
        ycoords = hole_wavefnc.infmass.ycoords;
        
        radius = 90*sqrt(3)/2*1.42;
        
        diameter_x = [mean(xcoords(1,:))-radius, mean(xcoords(1,:))+radius];
        diameter_y = [mean(ycoords(:,1)), mean(ycoords(:,1))];
        
        text_pos = [mean(xcoords(1,:))*0.76, mean(ycoords(:,1))*0.95];
        Level_list = [3.00735878351463e-34 4.46458986025152e-06 8.92917972050304e-06 1.33937695807546e-05 1.78583594410061e-05 2.23229493012576e-05 2.67875391615091e-05 3.12521290217606e-05 3.57167188820121e-05 4.01813087422637e-05 4.46458986025152e-05 4.91104884627667e-05 5.35750783230182e-05 5.80396681832697e-05 6.25042580435213e-05 6.69688479037728e-05 7.14334377640243e-05 7.58980276242758e-05 8.03626174845273e-05 8.48272073447788e-05 8.92917972050304e-05 9.37563870652819e-05 9.82209769255334e-05 0.000102685566785785 0.000107150156646036 0.000111614746506288];
    
        % plot the wave function
         contourf( xcoords(1:delta:end,1:delta:end), ycoords(1:delta:end,1:delta:end), abs(Wavefnc(1:delta:end,1:delta:end)).^2, ...
            'LevelList', Level_list, 'LineStyle', 'none',  'Parent', axes_geo1 );  
       
        % label the diameter
        plot( diameter_x, diameter_y, 'LineStyle', '-', 'LineWidth', 2, 'color', [0 0 0], 'Parent', axes_geo1 )
        text(text_pos(1), text_pos(2), [num2str(2*radius/10,'% 10.2f'),' nm'], 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo1);
        
        % creating energy text
        energy = hole_wavefnc.infmass.eigenval/homega;
        text(text_pos(1)-radius/5, text_pos(2)+2.1*radius, ['E = ',num2str(energy,'% 10.2f'),' $$\hbar\omega_c$$'], 'Interpreter', 'latex', 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo1);
        

        %************** wave function for smooth edge **********************
        
        axes_geo2 = axes( 'Parent',figure1, ...
                'CLim',[0 1e-4],...
                'Visible', 'off',...
                'FontSize', fontsize,...         
                'ylim', y_lim,...  
                'Units', 'Pixels', ...
                'FontName','Times New Roman');
        hold on; 
        axis equal;
              
        % selecting data to plot        
        Wavefnc = hole_wavefnc.smoothedge.Wavefnc{1};
        xcoords = hole_wavefnc.smoothedge.xcoords;
        ycoords = hole_wavefnc.smoothedge.ycoords;
        
        % plot the wave function
         contourf( xcoords(1:delta:end,1:delta:end), ycoords(1:delta:end,1:delta:end), abs(Wavefnc(1:delta:end,1:delta:end)).^2, ...
            'LevelList', Level_list, 'LineStyle', 'none',  'Parent', axes_geo2 ); 
       
       % label the diameter
       plot( diameter_x, diameter_y, 'LineStyle', '-', 'LineWidth', 2, 'color', [0 0 0], 'Parent', axes_geo2 )       
       text(text_pos(1), text_pos(2), [num2str(2*radius/10,'% 10.2f'),' nm'], 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo2);
       
       % creating energy text
       energy = hole_wavefnc.smoothedge.eigenval/homega;
       text(text_pos(1)-radius/5, text_pos(2)+2.1*radius, ['E = ',num2str(energy,'% 10.2f'),' $$\hbar\omega_c$$'], 'Interpreter', 'latex', 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo2);

        
       %************** wave function for rough edge **********************
        axes_geo3 = axes('Parent',figure1,...
                'CLim',[0 1e-4],...
                'Visible', 'off',...
                'FontSize', fontsize,...                      
                'ylim', y_lim,... 
                'Units', 'Pixels', ...
                'FontName','Times New Roman');
        hold on; 
        axis equal;
        
        
        Wavefnc = hole_wavefnc.roughedge.Wavefnc{1};
        xcoords = hole_wavefnc.roughedge.xcoords;
        ycoords = hole_wavefnc.roughedge.ycoords;
        
        % plot the wave function
         contourf( xcoords(1:delta:end,1:delta:end), ycoords(1:delta:end,1:delta:end), abs(Wavefnc(1:delta:end,1:delta:end)).^2, ...
           'LevelList', Level_list, 'LineStyle', 'none',  'Parent', axes_geo3 ); 
      
      % label the diameter
       plot( diameter_x, diameter_y, 'LineStyle', '-', 'LineWidth', 2, 'color', [0 0 0], 'Parent', axes_geo3 )     
       text(text_pos(1), text_pos(2), [num2str(2*radius/10,'% 10.2f'),' nm'], 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo3);
       
       % creating energy text
       energy = hole_wavefnc.roughedge.eigenval/homega;
       text(text_pos(1)-radius/5, text_pos(2)+2.1*radius, ['E = ',num2str(energy,'% 10.2f'),' $$\hbar\omega_c$$'], 'Interpreter', 'latex', 'FontSize', fontsize, 'FontName','Times New Roman',  'Parent', axes_geo3);
       
       % get the figure position
       figure_pos = get( figure1, 'Position' );
       
       
       cbar = colorbar('South', 'FontSize', fontsize, 'FontName','Times New Roman', 'XTickLabel', {'0', '',  'max'}, 'Units', 'Pixels');
       cbar_position = get(cbar, 'Position' );
       cbar_position(2) = 120;
       cbar_position(4) = 10;% = [30 10 figure_pos(3)/4.3 10];
       cbar_position(3) = figure_pos(3)/4.3;
       set( cbar, 'Position', cbar_position);
       
       
       
       % Set the position of the wave function with infinite mass
       Position_geo1 = get(axes_geo1, 'Position');
       Position_geo1(1) = 0;
       Position_geo1(2) = 0;
       Position_geo1(3) = figure_pos(3)/3;
       Position_geo1(4) = figure_pos(4);
       set(axes_geo1, 'Position', Position_geo1);    

       % Set the position of the wave function with smooth edge
       Position_geo2 = get(axes_geo2, 'Position');% - get(axes_geo2, 'TightInset') * [-1 0 1 0; 0 -1 0 1; 0 0 1 0; 0 0 0 1];
       Position_geo2(1) = figure_pos(3)/3;
       Position_geo2(2) = 0;
       Position_geo2(3) = figure_pos(3)/3;
       Position_geo2(4) = figure_pos(4);
       set(axes_geo2, 'Position', Position_geo2);  
       
       % Set the position of the wave function with rough edge
       Position_geo3 = get(axes_geo3, 'Position');
       Position_geo3(1) = 2*figure_pos(3)/3;
       Position_geo3(2) = 0;
       Position_geo3(3) = figure_pos(3)/3;
       Position_geo3(4) = figure_pos(4);
       set(axes_geo3, 'Position', Position_geo3);  
       
       
       print( '-depsc2', [outputdir,'/',outfilename,'.eps'] )       
        
       close( figure1 )
        
end  

%% setOutputDir
% @brief Sets the output directory
    function setOutputDir()
        resultsdir = 'results';
        mkdir(resultsdir );
        outputdir = resultsdir;        
    end




end