TMDC_Monolayer_Lead_Hamiltonians.m

Go to the documentation of this file.

%%   Eotvos Quantum Transport Utilities - TMDC_Monolayer_Lead_Hamiltonians
%    Copyright (C) 2009-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 lattices Lattices
%> @{
%> @file TMDC_Monolayer_Lead_Hamiltonians.m
%> @brief Class to create the Hamiltonian of one unit cell in a translational invariant lead made of TMDC_Monolayer lattice structure.
%> @} 
%> @brief Class to create the Hamiltonian of one unit cell in a translational invariant lead made of TMDC_Monolayer lattice structure.
%%
classdef TMDC_Monolayer_Lead_Hamiltonians


methods (Access = public)
%% TMDC_Monolayer_Hamiltonians
%> @brief Creates Hamiltonians H_0 and H_1 for TMDC_Monolayer lattice structure
%> @image html TMDC_monolayer.jpg
%> @image latex TMDC_monolayer.jpg
%> @param lead_param An instance of structure #Lattice_TMDC_Monolayer (or its subclass) containing the physical parameters.
%> @param M Number of sites in the cross section of the lead.
%> @return [1] The Hamiltonian of one slab in the ribbon. 
%> @return [2] The coupling between the slabs. 
%> @return [3] The transverse coupling between the slabs for transverse calculations. 
%> @return [4] A structure Coordinates containing the coordinates of the sites. 
    function [H0, H1, H1_transverse, H1_skew_left, cCoordinates] = TMDC_Monolayer_Hamiltonians( obj, lead_param, M )
        
        % check the structure containing the physical parameters
        supClasses = superclasses(lead_param);
        if sum( strcmp( supClasses, 'Lattice_TMDC_Monolayer' ) ) == 0
            error(['EQuUs:Lattices:', class(obj), ':TMDC_Monolayer_Hamiltonians'], 'Invalid type of the input parameter');   
        end
        
        if isempty(M)
            error(['EQuUs:Lattices:', class(obj), ':TMDC_Monolayer_Hamiltonians'], 'The input parameter M is empty')
        end        
            
        % defining orbital types according to Table II in <a href="https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.205108">PRB 92 205108</a>
        orbitals = {'d0_xz', 'd0_yz', 'p0_z', 'p0_x', 'p0_y', 'de_z2', 'de_xy', 'de_x2y2', 'pe_z', 'pe_x', 'pe_y'};
        
        % number of orbitals
        orbitals_num = 11;
        
        % coordinates shift of the atomic sites 
        a1 = [1; 0];
        a2 = [-0.5; sqrt(3)/2 ];
        coordinates_shift = (2*a1+a2)/3;
        coordinates_shifts = { [0;0], [0;0], coordinates_shift, coordinates_shift, coordinates_shift, [0;0], [0;0], [0;0], coordinates_shift, coordinates_shift, coordinates_shift };
        
        % creating cell arrays for the Hamiltonians
        H0_orbitals = cell( orbitals_num, 1);
        H1_orbitals = cell( orbitals_num, 1);
        H1_orbitals_transverse = cell( orbitals_num, 1);
        H1_orbitals_skew_left = cell( orbitals_num, 1);
        
        % creating cell array of the coordinates
        cCoordinates_orbitals = cell( orbitals_num, 1);
        
        % obtaining derived hopping parameters t_2__i_i
        t_2__hoppings = lead_param.Calc_t_2__i_i();
        
        % creating triangular lattice for each atomic orbital
        for idx = 1:orbitals_num
    
            parameter_extension = ['__', num2str(idx), '_', num2str(idx)];
            epsilon = lead_param.(['epsilon', num2str(idx)]);
            t_1 = lead_param.(['t_1', parameter_extension]);
            t_2 = t_2__hoppings.(['t_2', parameter_extension]);
            [H0_orbitals{idx}, H1_orbitals{idx}, H1_orbitals_transverse{idx}, H1_orbitals_skew_left{idx}, cCoordinates_orbitals{idx}] = obj.Triangle_Lattice_Hamiltonians( epsilon, t_1, t_2, M );
            cCoordinates_orbitals{idx}.OrbitalTypes = cell( length(cCoordinates_orbitals{idx}.x), 1);
            cCoordinates_orbitals{idx}.OrbitalTypes(:) = {orbitals{idx}};
            
            % shift the given sublattice coordinates
            cCoordinates_orbitals{idx} = cCoordinates_orbitals{idx}.Shift( coordinates_shifts{idx} );
            
        end        
        
        % creating the big Hamiltonian from the block-diagonal parts
        H0 = blkdiag(H0_orbitals{:});       
        H1 = blkdiag(H1_orbitals{:});
        H1_transverse = blkdiag(H1_orbitals_transverse{:});
        H1_skew_left = blkdiag(H1_orbitals_skew_left{:});
        
        % combining coordinates structures
        cCoordinates = cCoordinates_orbitals{1};
        for idx = 2:orbitals_num
            cCoordinates = cCoordinates.Combine( cCoordinates_orbitals{idx} );
        end
        
        % setting the lattice constant
        cCoordinates.LatticeConstant = lead_param.a;
        
        
        
        
        % ********** adding off-diagonal coupling elements between the same atom types *********
        
        % obtaining derived hopping parameters t_2__i_j
        t_2__hoppings = lead_param.Calc_t_2__i_j();
        
        % obtaining derived hopping parameters t_3__i_j
        t_3__hoppings = lead_param.Calc_t_3__i_j();
        
        % M_2 symmetry (+), EQ (5) in <a href="https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.205108">PRB 92 205108</a>
        orbital_indexes1 = [3, 6, 9];
        orbital_indexes2 = [5, 8, 11]; 
        for idx = 1:length( orbital_indexes1 )
            orbital_index1 = orbital_indexes1(idx);
            orbital_index2 = orbital_indexes2(idx);
            row_indexes =  (orbital_index1-1)*M + (1:M);
            col_indexes =  (orbital_index2-1)*M + (1:M);
            
            % +delta_1 direction
            hopping_ampliude_4 = lead_param.(['t_1__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( col_indexes, row_indexes ) = H1( col_indexes, row_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            
            % +delta_2 direction
            hopping_ampliude_4 = t_2__hoppings.(['t_2__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( col_indexes, row_indexes ) = H1( col_indexes, row_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_4, M, M );
            
            hopping_ampliude_3 = t_3__hoppings.(['t_3__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_3, M, M );
            
            H1_skew_left( col_indexes, row_indexes ) = H1_skew_left( col_indexes, row_indexes ) + sparse( M, 1, hopping_ampliude_4, M, M);
            H1_skew_left( row_indexes, col_indexes ) = H1_skew_left( row_indexes, col_indexes ) + sparse( M, 1, hopping_ampliude_3, M, M);
            
            % +delta_3 directions
            hopping_ampliude_4 = t_2__hoppings.(['t_2__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            hopping_ampliude_3 = t_3__hoppings.(['t_3__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_4, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_3, M, M );
            H1_transverse( col_indexes, row_indexes ) =  H1_transverse( col_indexes, row_indexes ) + sparse( M, 1, hopping_ampliude_4, M, M);
            H1_transverse( row_indexes, col_indexes ) =  H1_transverse( row_indexes, col_indexes ) + sparse( M, 1, hopping_ampliude_3, M, M);
            
            % - delta_3 direction
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) + sparse( 2:M, 1:M-1, hopping_ampliude_3, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) + sparse( 2:M, 1:M-1, hopping_ampliude_4, M, M );  
        end
        
        
        
        
        
        % M_2 symmetry (-), EQ (6) in <a href="https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.205108">PRB 92 205108</a>
        orbital_indexes1 = [1, 3, 4, 6, 7, 9, 10];
        orbital_indexes2 = [2, 4, 5, 7, 8, 10, 11];        
        for idx = 1:length( orbital_indexes1 )
            orbital_index1 = orbital_indexes1(idx);
            orbital_index2 = orbital_indexes2(idx);
            row_indexes =  (orbital_index1-1)*M + (1:M);
            col_indexes =  (orbital_index2-1)*M + (1:M);
            
            % -delta_1 direction
            hopping_ampliude_1 = lead_param.(['t_1__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( col_indexes, row_indexes ) = H1( col_indexes, row_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_1, M, M );
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) - sparse( 1:M, 1:M, hopping_ampliude_1, M, M );
            
            % +delta_2 direction
            hopping_ampliude_2 = t_2__hoppings.(['t_2__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( col_indexes, row_indexes ) = H1( col_indexes, row_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_2, M, M );
            
            hopping_ampliude_3 = t_3__hoppings.(['t_3__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) - sparse( 1:M-1, 2:M, hopping_ampliude_3, M, M );
            
            H1_skew_left( col_indexes, row_indexes ) = H1_skew_left( col_indexes, row_indexes ) + sparse( M, 1, hopping_ampliude_2, M, M);
            H1_skew_left( row_indexes, col_indexes ) = H1_skew_left( row_indexes, col_indexes ) - sparse( M, 1, hopping_ampliude_3, M, M);
            
            % +delta_3 directions
            hopping_ampliude_2 = t_2__hoppings.(['t_2__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            hopping_ampliude_3 = t_3__hoppings.(['t_3__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) - sparse( 1:M-1, 2:M, hopping_ampliude_2, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_3, M, M );
            H1_transverse( col_indexes, row_indexes ) =  H1_transverse( col_indexes, row_indexes ) - sparse( M, 1, hopping_ampliude_2, M, M);
            H1_transverse( row_indexes, col_indexes ) =  H1_transverse( row_indexes, col_indexes ) + sparse( M, 1, hopping_ampliude_3, M, M);
            
            % - delta_3 direction
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) + sparse( 2:M, 1:M-1, hopping_ampliude_3, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) - sparse( 2:M, 1:M-1, hopping_ampliude_2, M, M );
        end
        
        
        
        % ********** adding off-diagonal coupling elements between the different atom types *********
        
        % obtaining derived hopping parameters t_4__i_j
        t_4__hoppings = lead_param.Calc_t_4__i_j();
        
        % M_2 symmetry (+), EQ (7) in <a href="https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.205108">PRB 92 205108</a>
        orbital_indexes1 = [3, 5, 4, 10, 9, 11, 10];
        orbital_indexes2 = [1, 1, 2, 6, 7, 7, 8]; 
        for idx = 1:length( orbital_indexes1 )
            orbital_index1 = orbital_indexes1(idx);
            orbital_index2 = orbital_indexes2(idx);
            row_indexes =  (orbital_index1-1)*M + (1:M);
            col_indexes =  (orbital_index2-1)*M + (1:M);
            
            % +/-delta_4 direction
            hopping_ampliude_4 = t_4__hoppings.(['t_4__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            
            % + delta_6 direction
            hopping_ampliude_4 = t_4__hoppings.(['t_4__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) - sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            
        end
        
        
        
        
        % M_2 symmetry (-), EQ (8) in <a href="https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.205108">PRB 92 205108</a>
        orbital_indexes1 = [4, 3, 5, 9, 11, 10, 9, 11];
        orbital_indexes2 = [1, 2, 2, 6, 6, 7, 8, 8]; 
        for idx = 1:length( orbital_indexes1 )
            orbital_index1 = orbital_indexes1(idx);
            orbital_index2 = orbital_indexes2(idx);
            row_indexes =  (orbital_index1-1)*M + (1:M);
            col_indexes =  (orbital_index2-1)*M + (1:M);
            
            % +/-delta_4 direction
            hopping_ampliude_4 = t_4__hoppings.(['t_4__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H0( col_indexes, row_indexes ) = H0( col_indexes, row_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            H0( row_indexes, col_indexes ) = H0( row_indexes, col_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            
            % + delta_6 direction
            hopping_ampliude_4 = t_4__hoppings.(['t_4__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) + sparse( 1:M, 1:M, hopping_ampliude_4, M, M );
            
            % + delta_5 direction
            hopping_ampliude_5 = lead_param.(['t_5__', num2str(orbital_index1), '_', num2str(orbital_index2)]);
            H1( row_indexes, col_indexes ) = H1( row_indexes, col_indexes ) + sparse( 1:M-1, 2:M, hopping_ampliude_5, M, M );
            H1_skew_left( row_indexes, col_indexes ) =  H1_skew_left( row_indexes, col_indexes ) + sparse( M, 1, hopping_ampliude_5, M, M);
            
        end
        
    
    end
    
    
end % end public methods


methods (Access = protected, Static = true)
    

%% Triangle_Lattice_Hamiltonians
%> @brief Creates Hamiltonians H_0 and H_1 for a simple triangle lattice
%> @param epsilon The on-site energy
%> @param t_1 ?????????
%> @param t_2 ?????????
%> @param M Number of sites in the cross section of the lead.
%> @return [1] The Hamiltonian of one slab in the ribbon. 
%> @return [2] The coupling between the slabs. 
%> @return [3] The transverse coupling between the slabs for transverse calculations. 
%> @return [4] A structure Coordinates containing the coordinates of the sites. 
    function [H0, H1, H1_transverse, H1_skew_left, cCoordinates] = Triangle_Lattice_Hamiltonians( epsilon, t_1, t_2, M )
                
        if isempty(M)
            error(['EQuUs:TMDC_Monolayer_Lead_Hamiltonians:Triangle_Lattice_Hamiltonians'], 'The input parameter M is empty')
        end      
    
        H0 = sparse( 1:M, 1:M, epsilon, M, M) + ... on-site terms
            sparse( 2:M, 1:M-1, t_2, M, M) + ... hopping term (-delta_3 direction)
            sparse( 1:M-1, 2:M, t_2, M, M); % hopping term (delta_3 direction)
            
        H1 = sparse( 1:M, 1:M, t_1, M, M) + ... hopping terms  (delta_1 direction)
            sparse( 1:M-1, 2:M, t_2, M, M); % hopping terms  (delta_2 direction)
            
        % creating coordinates
        cCoordinates = structures('coordinates');
        cCoordinates.LatticeConstant = 1;
        cCoordinates.a = [1; 0];
        cCoordinates.b = [-0.5; sqrt(3)/2 ];
        
        cCoordinates.x = cCoordinates.b(1)*(1:M)';
        cCoordinates.y = cCoordinates.b(2)*(1:M)';
            
        % Transverse coupling    
        H1_transverse = sparse( M, 1, t_2, M, M);  % (delta_3 direction)
        
        % skew coupling in the left direction
        H1_skew_left = sparse( M, 1, t_2, M, M); % (delta_2 direction)
    
    end
    
    
end % end static methods

end