2 % Copyright (C) 2016 Peter Rakyta, Ph.D.
4 % This program is free software: you can redistribute it and/or modify
5 % it under the terms of the GNU General Public License as published by
6 % the Free Software Foundation, either version 3 of the License, or
7 % (at your option) any later version.
9 % This program is distributed in the hope that it will be useful,
10 % but WITHOUT ANY WARRANTY; without even the implied warranty of
11 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 % GNU General Public License
for more details.
14 % You should have received a copy of the GNU General Public License
15 % along with
this program. If not, see http:
17 %> @addtogroup basic Basic Functionalities
20 %> @brief Class to create and store Hamiltonian of the translational invariant leads.
22 %> @brief Class to create and store Hamiltonian of the translational invariant leads.
24 %> EQuUs v4.8 or later
29 properties ( Access =
protected )
30 %> An instance of the structure
param 32 %> The orientation of the lead. Set +1 is the
"incoming" direction of the propagating states is defined in the +x or +y direction, and
"-1" otherwise.
34 %> The
id number of the current lead.
36 %> Function
handle to calculate the infinite Greens
function.
38 %> The number of the
sites in the cross section.
42 %> The tranverse momentum
for transverse computations.
44 %> List of
sites in the unit cell that should be kept after decimation.
46 %> An instance of the structure coordinates.
48 %> K0=H0-E*S0, see Eq (4) of PRB 78, 035407
50 %> K1=H1-E*S1, see Eq (4) of PRB 78, 035407
52 %> K1adj=H1adj-E*S1', see Eq (4) of PRB 78, 035407
54 %> K1_transverse=H1_transverse-E*S1_transverse
56 %> The Hamiltonian of a unit cell.
58 %> The coupling Hamiltonian between the unit cells.
60 %> The coupling Hamiltonian between the unit cells in the opposite direction as H1. (For complex energies they differ from each other.)
64 %> The transverse coupling between the slabs for transverse calculations.
66 %> The overlap integrals of a unit cell.
68 %> The overlap integrals between the unit cells.
70 %> The adjungate of the overlap integrals between the unit cells.
72 %> The overlap integrals between the slabs for transverse calculations.
74 %> The matrix of the
Peierls phases in the unit cell.
76 %> The matrix of the
Peierls phases in the coupling matrix between the unit cells.
78 %> The matrix of the
Peierls phases in the transverse coupling matrix between the unit cells.
80 %> A logical value. True if the
Hamiltonians were created, false otherwise.
82 %> A logical value. True if the
Hamiltonians were decimated, false otherwise.
84 %> A logical value. True if the overlap integrals were applied, false otherwise.
86 %> A logical value. True if magnetic field was applied in the
Hamiltonians, false otherwise.
88 %> A logical value. True if a gauge transformation was incorporated into the Hamiltonians or false otherwise.
89 GaugeTransformationApplied
90 %> list of optional parameters (see http:
95 methods ( Access = public )
96 %% constructorof the class
97 %> @brief Constructor of the class.
100 %> @
param varargin Cell array of optional parameters. See
#InputParsing for details. 101 %> @
return An instance of the
class 105 obj.varargin = varargin;
114 %% ApplyOverlapMatrices
115 %> @brief Applies the overlap matrices to the Hamiltonians: K = H-ES
116 %> @param E The energy value.
117 function ApplyOverlapMatrices(obj, E)
119 if obj.OverlapApplied
124 if ~isempty( obj.S0 )
125 obj.K0 = obj.H0 - E*obj.S0;
128 obj.K0 = obj.H0 - speye(size(obj.H0))*E;
131 if ~isempty( obj.S1 )
132 obj.K1 = obj.H1 - E*obj.S1;
133 obj.K1adj = obj.H1adj - E*obj.S1';
136 obj.K1adj = obj.H1adj;
140 if ~isempty( obj.S1_transverse )
141 obj.K1_transverse = obj.H1_transverse - E*obj.S1_transverse;
142 elseif ~isempty( obj.H1_transverse )
143 obj.K1_transverse = obj.H1_transverse;
146 obj.OverlapApplied = true;
151 %> @brief Creates the Hamiltonians H_0 and H_1 of the lead. The created Hamiltonians are stored by within the
object.
152 %> @param varargin Cell array of optional parameters (https:
153 %> @param 'toSave' Logical value. If true, the created Hamiltonians are saved into a file 'Hamiltoni_Lead_' + num2str(Hanyadik_Lead) + '.mat'.
154 %> @param 'CustomHamiltonian' An instance of class
#Custom_Hamiltonians describing external source of Hamiltonians. 158 p.addParameter(
'toSave', 0);
159 p.addParameter(
'CustomHamiltonian', []);
160 p.parse(varargin{:});
162 toSave = p.Results.toSave;
163 CustomHamiltonian = p.Results.CustomHamiltonian;
168 if strcmpi(obj.Opt.Lattice_Type,
'Square')
170 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.SquareLattice_Lead_Hamiltonians(obj.params, obj.M,
'q', obj.q);
172 obj.kulso_szabfokok = 1:obj.M; 173 elseif strcmpi(obj.Opt.Lattice_Type, 'SSH
') 174 createH = Square_Lead_Hamiltonians(); 175 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.SSH_Lead_Hamiltonians(obj.params, 'q
', obj.q); 177 obj.kulso_szabfokok = 1;
178 elseif strcmp(obj.Opt.Lattice_Type,
'Lieb')
180 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Lieb_Lead_Hamiltonians(obj.params, obj.M,
'q', obj.q);
182 obj.M = size(obj.H0,1); 183 obj.kulso_szabfokok = [];%1:3:3*obj.M-2; 184 elseif strcmp(obj.Opt.Lattice_Type, 'BiTeI
') 185 createH = BiTeI_Lead_Hamiltonians(); 186 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.BiTeILattice_Lead_Hamiltonians(obj.params, obj.M, 'q
', obj.q); 188 obj.M = size(obj.H0,1);
189 obj.kulso_szabfokok = [];%sort([1:4:2*obj.M-1, 2:4:2*obj.M],
'ascend');%[];
190 elseif strcmp(obj.Opt.Lattice_Type,
'Graphene') || strcmpi(obj.Opt.Lattice_Type,
'H')
192 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Graphene_Lead_Hamiltonians(obj.params, obj.M, obj.End_Type,
'q', obj.q);
194 obj.kulso_szabfokok = (1:obj.M); 195 elseif strcmpi(obj.Opt.Lattice_Type, 'Graphene_Bilayer
') 196 createH = Hex_Lead_Hamiltonians(); 197 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Graphene_Bilayer_Lead_Hamiltonians(obj.params, obj.M, obj.End_Type, 'q
', obj.q); 199 obj.kulso_szabfokok = [1:obj.M, size(obj.H0,1)/2 + (1:obj.M)];
201 elseif strcmpi(obj.Opt.Lattice_Type,
'Graphene_Bilayer_2')
203 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Graphene_Bilayer_Lead_Hamiltonians2(obj.params, obj.M, obj.End_Type,
'q', obj.q);
205 obj.kulso_szabfokok = [1:obj.M, size(obj.H0,1)/2 + (1:obj.M)]; 207 elseif strcmpi(obj.Opt.Lattice_Type, 'Graphene_Bilayer_3
') 208 createH = Hex_Lead_Hamiltonians(); 209 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Graphene_Bilayer_Lead_Hamiltonians3(obj.params, obj.M, obj.End_Type, 'q
', obj.q); 211 obj.kulso_szabfokok = [1:obj.M+1, size(obj.H0,1)/2 + (1:obj.M+1)];
213 elseif strcmpi(obj.Opt.Lattice_Type,
'Silicene')
215 [obj.H0,obj.H1,obj.H1_transverse,obj.coordinates] = createH.Silicene_Lead_Hamiltonians(obj.params, obj.M, obj.End_Type,
'q', obj.q);
217 obj.kulso_szabfokok = [1:obj.M, size(obj.H0,1)/2 + (1:obj.M)]; 219 elseif ~isempty( obj.Opt.custom_Hamiltonians ) 220 if isempty( CustomHamiltonian ) 221 CustomHamiltonian = Custom_Hamiltonians( obj.Opt, obj.param ); 224 if ~CustomHamiltonian.Read( 'Hamiltonians_loaded
' ) 225 CustomHamiltonian.LoadHamiltonians(); 228 obj.coordinates = CustomHamiltonian.Read( 'coordinates
' ); 229 obj.coordinates = obj.coordinates{obj.Hanyadik_Lead}; 230 obj.H0 = CustomHamiltonian.Read( 'H0
' ); 231 obj.H0 = obj.H0{obj.Hanyadik_Lead}; 232 obj.H1 = CustomHamiltonian.Read( 'H1
' ); 233 obj.H1 = obj.H1{obj.Hanyadik_Lead}; 235 obj.H1_transverse = CustomHamiltonian.Read(
'H1_transverse' );
236 obj.H1_transverse = obj.H1_transverse{obj.Hanyadik_Lead};
237 obj.S0 = CustomHamiltonian.Read(
'S0' );
239 obj.S0 = obj.S0{obj.Hanyadik_Lead};
241 obj.S1 = CustomHamiltonian.Read(
'S1' );
243 obj.S1 = obj.S1{obj.Hanyadik_Lead};
246 obj.S1_transverse = CustomHamiltonian.Read( 'S1_transverse
' ); 247 if iscell( obj.S1_transverse ) 248 obj.S1_transverse = obj.S1_transverse{obj.Hanyadik_Lead}; 250 obj.M = size(obj.H0,1); 251 %obj.kulso_szabfokok = 1:obj.M; 253 error(['EQuUs:
', class(obj), ':
CreateHamiltonians'], 'Unrecognized lattice type, or a valid custom source
for the Hamiltonians was not
set.
') 256 obj.Transform2Spin(); 263 obj.HamiltoniansCreated = true; 264 obj.OverlapApplied = false; 265 obj.HamiltoniansDecimated = false; 266 obj.MagneticFieldApplied = false; 267 obj.GaugeTransformationApplied = false; 273 %> @brief Transforms the Hamiltonians and the overlap matrices to include electron spin. 274 function Transform2Spin( obj ) 276 if isempty(obj.Opt.Spin) || ~obj.Opt.Spin 280 if ~isempty( obj.coordinates.spinup ) 281 % spin is also included in the model 285 fnames = fieldnames( obj.coordinates ); 286 for idx = 1:length(fnames) 288 if strcmp(fname, 'a
') || strcmp(fname, 'b
') 291 obj.coordinates.(fname) = [ obj.coordinates.(fname); obj.coordinates.(fname) ]; 294 obj.coordinates.spinup = [ true(size(obj.H0,1),1); false(size(obj.H0,1),1)]; 296 obj.kulso_szabfokok = [obj.kulso_szabfokok, obj.kulso_szabfokok+size(obj.H0,1)]; 299 % transforming the Hamiltonians 300 obj.H0 = [obj.H0, sparse([],[],[], size(obj.H0,1), size(obj.H0,2)); sparse([],[],[], size(obj.H0,1), size(obj.H0,2)), obj.H0]; 301 obj.H1 = [obj.H1, sparse([],[],[], size(obj.H1,1), size(obj.H1,2)); sparse([],[],[], size(obj.H1,1), size(obj.H1,2)), obj.H1]; 302 obj.H1adj = [obj.H1adj, sparse([],[],[], size(obj.H1adj,1), size(obj.H1adj,2)); sparse([],[],[], size(obj.H1adj,1), size(obj.H1adj,2)), obj.H1adj]; 304 obj.H1_transverse = [obj.H1_transverse, sparse([],[],[], size(obj.H1_transverse,1), size(obj.H1_transverse,2)); ... 305 sparse([],[],[], size(obj.H1_transverse,1), size(obj.H1_transverse,2)), obj.H1_transverse]; 307 % transforming th eoverlap integrals 308 if ~isempty( obj.S0 ) 309 obj.S0 = [obj.S0, sparse([], [], [], size(obj.S0,1), size(obj.S0,2)); ... 310 sparse([], [], [], size(obj.S0,1), size(obj.S0,2)), obj.S0]; 314 obj.S1 = [obj.S1, sparse([], [], [], size(obj.S1,1), size(obj.S1,2)); ... 315 sparse([], [], [], size(obj.S1,1), size(obj.S1,2)), obj.S1]; 318 if ~isempty( obj.S1_transverse ) 319 obj.S1_transverse = [obj.S1_transverse, sparse([], [], [], size(obj.S1_transverse,1), size(obj.S1_transverse,2)); ... 320 sparse([], [], [], size(obj.S1_transverse,2), size(obj.S1_transverse,1)), obj.S1_transverse]; 329 %> @brief Transforms the Hamiltonians and the overlap matrices into the BdG model in the Nambu space representation according to 330 %> <a href="http://iopscience.iop.org/article/10.1088/1367-2630/9/8/278/meta">New Journal of Physics 9 (2007) 278</a>. 331 %> It is assumed, that the Hamiltonian is already transfromed to the grand canonical operator: \f$ \hat{H} \rightarrow \hat{H} - E_F\hat{N}\f$ 332 function Transform2BdG( obj ) 334 if isempty(obj.Opt.BdG) || ~obj.Opt.BdG 335 obj.display(['EQuUs:
', class(obj), ':Transform2BdG: BdG option is not
set to
true in the computational parameters.
']) 339 if ~isempty( obj.coordinates.BdG_u ) 340 % already transformed into BdG model 341 obj.display(['EQuUs:
', class(obj), ':Transform2BdG: Hamiltonians already transformed intt the BdG model.
']) 345 fnames = fieldnames( obj.coordinates ); 346 for idx = 1:length(fnames) 348 if strcmp(fname, 'a
') || strcmp(fname, 'b
') 351 obj.coordinates.(fname) = [ obj.coordinates.(fname); obj.coordinates.(fname) ]; 354 obj.coordinates.BdG_u = [ true(size(obj.H0,1),1); false(size(obj.H0,1),1)]; 356 obj.kulso_szabfokok = [obj.kulso_szabfokok, obj.kulso_szabfokok+size(obj.H0,1)]; 358 pair_potential = obj.params.pair_potential; 360 % transforming the Hamiltonians 362 S0 = speye(size(obj.H0)); 367 obj.H0 = [obj.H0, S0*pair_potential; S0*conj(pair_potential), -conj(obj.H0)]; 370 S1 = sparse([],[],[], size(obj.H1,1), size(obj.H1,2)); 375 obj.H1 = [obj.H1, S1*pair_potential; S1*conj(pair_potential), -conj(obj.H1)]; 376 obj.H1adj = [obj.H1adj, S1'*pair_potential; S1
'*conj(pair_potential), -conj(obj.H1adj)]; 378 if isempty(obj.S1_transverse) 379 S1_transverse = sparse([],[],[], size(obj.H1_transverse,1), size(obj.H1_transverse,2)); 381 S1_transverse = obj.S1_transverse; 384 obj.H1_transverse = [obj.H1_transverse, S1_transverse*pair_potential; S1_transverse*conj(pair_potential), -conj(obj.H1_transverse)]; 386 % transforming the overlap integrals 387 if ~isempty( obj.S0 ) 388 obj.S0 = [obj.S0, sparse([], [], [], size(obj.S0,1), size(obj.S0,2)); sparse([], [], [], size(obj.S0,1), size(obj.S0,2)), obj.S0]; 392 obj.S1 = [obj.S1, sparse([], [], [], size(obj.S1,1), size(obj.S1,2)); sparse([], [], [], size(obj.S1,1), size(obj.S1,2)), obj.S1]; 395 if ~isempty( obj.S1_transverse ) 396 obj.S1_transverse = [obj.S1_transverse, sparse([], [], [], size(obj.S1_transverse,1), size(obj.S1_transverse,2)); ... 397 sparse([], [], [], size(obj.S1_transverse,2), size(obj.S1_transverse,1)), obj.S1_transverse]; 400 % the number of sites in the cross section becomes twice as many as in the normal case 406 %> @brief Calculates the band structure of the lead. 407 %> @param varargin Cell array of optional parameters: 408 %> @param 'toPlot
' Set 1 in order to plot the calculated spectrum, 0 (default) otherwise 409 %> @param 'ka_min
' The lower bound of the wave numbers. (Default is -pi.) 410 %> @param 'ka_max
' The upper bound of the wave numbers. (Default is pi.) 411 %> @param 'ka_num
' The number of wave number points involved in the calculations. (Default is 300.) 412 %> @param 'ka_vec
' One dimensional array of the k-pints. (Overrides previous attributes related to the k-vector array.) 413 %> @param 'center
' The calculated energy eigenvalues are centered around this value. (Default is 0.001.) 414 %> @param 'db
' The number of the calculated eigenvalues. 415 %> @param 'offset
' Offset value to shift the spectrum along the energy axis. 416 %> @param 'calcWaveFnc
' Logical value. Set true to calculate also the wave functions, or false (default) otherwise. 417 %> @return [1] ka_num x 2 array of the calculated spactrum. In the first column are the k-points, whil ein the second columns are the calculated energy points. 418 %> @return [2] The calculated wave functions stored in a structure #WaveFnc. 419 function [spectrum, WaveFnc] = CalcSpektrum( obj, varargin ) 422 p.addParameter('toPlot
', 0, @isscalar); 423 p.addParameter('ka_min
', -pi, @isscalar); 424 p.addParameter('ka_max
', pi, @isscalar); 425 p.addParameter('ka_num
', 300, @isscalar); 426 p.addParameter('ka_vec
', [] ); 427 p.addParameter('center
', 0.001); 428 p.addParameter('db
', min([10,size(obj.H0,1)]), @isscalar); 429 p.addParameter('offset
', mean(diag(obj.H0)) ); %Offset value to shift the spectrum along the energy axis 430 p.addParameter('calcWaveFnc
', false ); %Logical value. Set true to calculate also the wave functions, or false (default) otherwise. 432 p.parse(varargin{:}); 433 toPlot = p.Results.toPlot; 434 ka_min = p.Results.ka_min; 435 ka_max = p.Results.ka_max; 436 ka_num = p.Results.ka_num; 437 ka_vec = p.Results.ka_vec; 438 center = p.Results.center; 440 offset = p.Results.offset; 441 calcWaveFnc = p.Results.calcWaveFnc; 444 % check wheter Hamiltonians are decimated or not 445 if obj.HamiltoniansDecimated 446 obj.display(['EQuUs:
', class(obj), ':CalcSpektrum: Hamiltonians are decimated. Please recreate Hamiltonians to calculate spectrum.
'], 1) 451 % check the number of eigenvalues 452 db = min([ db, size(obj.H0,1)]); 455 obj.display('Calculating spectrum
') 457 % discrete increment of the wavenumber array 458 deltak = (ka_max-ka_min)/ka_num; 460 % the startpoint of the wavenumbers 463 % allocating temporary matrices 465 H0loc = obj.H0 - eye(size(obj.H0))*offset; 466 H1_transverse_loc = obj.H1_transverse; 469 S1_transverse_loc = obj.S1_transverse; 471 % the transverse momentum number 475 % creating the one-dimensional array for the wave numbers if not given 477 ka_vec = ka_min:deltak:ka_max; 480 % allocating arrays for the results 481 spectrum = cell(length(ka_vec),1); 483 WaveFnc = structures('WaveFnc'); 488 % calculating the spectrum 489 for idx=1:length(ka_vec) 490 % obtaining the k-dependent effective Hamiltonian 491 H = secular_H( H0loc, H1loc, H1_transverse_loc, ka_vec(idx)); 493 if isempty( S0_loc ) && isempty( S1_loc ) 495 % calculations including the wavefunctions 497 [WaveFnc_tmp, E] = eigs(H, db, center); 499 [WaveFnc_tmp, E] = eig(H); 502 WaveFnc(idx).WaveFnc = WaveFnc_tmp; 504 WaveFnc(idx).ka = ka_vec(idx); 506 % only eigenvalues are calculated 508 E = eigs(H, db, center); 514 % calculations including the overlap matrices 515 S = secular_H( S0_loc, S1_loc, S1_transverse_loc, ka_vec(idx)); 517 % calculations including the wavefunctions 519 [WaveFnc_tmp, E] = eigs(H, S, db, center); 521 [WaveFnc_tmp, E] = eig(H, S); 524 WaveFnc.WaveFnc{idx} = WaveFnc_tmp; 526 WaveFnc.ka(idx) = ka_vec(idx); 528 % only eigenvalues are calculated 530 E = eigs(H, S, db, center); 537 spectrum{idx} = [ones(length(E),1)*ka_vec(idx), E]; 540 % reorganize the calculated data 541 spectrum = cell2mat(spectrum); 544 obj.display('Spectrum calculated
') 546 % check whether to plot the spectrum 551 % plot the calculated spectrum 555 spectrum(:,2) = real(spectrum(:,2)); 556 x_lim = [min(spectrum(:,1)) max(spectrum(:,1))]; 557 y_lim = [min(spectrum(:,2)) max(spectrum(:,2))]; 559 Position = [0.5 0.58 0.33 0.4]; 560 axes_spectrum = axes('Parent
',figure1, 'Position
', Position,... 562 'FontSize
', fontsize,... 564 'ylim
', y_lim,... 'XTick
', XTick,... 'YTick
', YTick,... 566 'FontName
','Times New Roman
'); 569 plot(spectrum(:,1), spectrum(:,2),'.
', 'MarkerSize
', 4, 'Parent
', axes_spectrum ) 572 xlabel_position = [0 0 0]; 573 xlabel('ka
','FontSize
', fontsize,'FontName
','Times New Roman
', 'Parent
', axes_spectrum); 574 xlabel_handle = get(axes_spectrum,'XLabel
'); 575 set(xlabel_handle, 'Position
', get(xlabel_handle, 'Position
') + xlabel_position); 578 ylabel_position = [0 0 0]; 579 ylabel('E [eV]
','FontSize
', fontsize,'FontName
','Times New Roman
', 'Parent
', axes_spectrum); 580 ylabel_handle = get(axes_spectrum,'YLabel
'); 581 set(ylabel_handle, 'Position
', get(ylabel_handle, 'Position
') + ylabel_position); 584 % -------------------------------------- 586 % Hamiltonian for the secular equation 587 function H = secular_H(H0,H1, H1_transverse, k) 589 if ~isempty(q) && ~obj.HamiltoniansDecimated 590 H0 = H0 + H1_transverse*diag(exp(-1i*q)) + H1_transverse'*diag(exp(1i*q));
592 H = H0 + H1*exp(1i*k) + H1
'*exp(-1i*k); 596 % end nested functions 602 %> @brief Save Lead Hamiltonians into a file 'Hamiltoni_Lead_
' + num2str(Hanyadik_Lead) + '.mat
'. 603 function saveLeads( obj ) 604 save(['Hamiltoni_Lead_
',num2str(obj.Hanyadik_Lead),'.mat
'], 'H0
', 'H1
', 'kulso_szabfokok
', 'Lead_Orientation
', 'fazis_mtx_H0
', 'fazis_mtx_H1
', 'params
'); 608 %> @brief Shifts the coordinates of the sites by an integer multiple of the lattice vector #Coordinates.a. 609 %> @param shift Integer by which the coordinates are shifted. 610 function ShiftCoordinates( obj, shift ) 611 if isempty(shift) || shift == 0 615 % shifting the coordinates along the translational invariant direction 616 obj.coordinates = obj.coordinates.Shift( shift*obj.coordinates.a ); 623 %> @brief Shifts the on-site energies in the leads by a given energy. 624 %> @param Energy The enrgy value. 625 function ShiftLead( obj, Energy ) 626 obj.H0 = obj.H0 + sparse(1:size(obj.H0,1), 1:size(obj.H0,2),Energy,size(obj.H0,1),size(obj.H0,2)); 627 obj.params.epsilon = obj.params.epsilon + Energy; 630 obj.OverlapApplied = false; 634 %> @brief Adds on-site potential to the Hamiltonian H0. 635 %> @param V The potential calculated on the sites. 636 function AddPotential( obj, V ) 638 if ( size(V,2) == 1 ) || ( size(V,1) == 1 ) 639 if length(V) == size(obj.H0,1) 640 obj.H0 = obj.H0 + sparse(1:size(obj.H0,1), 1:size(obj.H0,2),V,size(obj.H0,1),size(obj.H0,2)); 642 obj.OverlapApplied = false; 644 disp(' Wrong dimension of potential: 1
') 647 elseif norm(size(V) - size(obj.H0)) < 1e-6 650 obj.OverlapApplied = false; 652 disp(' Wrong dimension of potential: 2
') 656 obj.display(' Potential added to lead
') 660 %> @brief Test, whether the lead is in the superconducting phase or not. 661 %> @return True if the lead is superconducting, false otherwise. 662 function ret = isSuperconducting( obj ) 669 if abs(obj.params.pair_potential) >= 1e-10 681 %> @brief Creates a clone of the present class. 682 %> @return Returns with the cloned object. 683 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html): 684 %> @param 'empty
' Set true to create an empty clone, or false (default) to clone all atributes. 685 function ret = CreateClone( obj, varargin ) 688 p.addParameter('empty
', false ); %Logical value. Set true for creating an empty class, or false (default) otherwise. 690 p.parse(varargin{:}); 691 empty = p.Results.empty; 693 ret = CreateLeadHamiltonians(obj.Opt, obj.param,... 694 'Hanyadik_Lead
', obj.Hanyadik_Lead,... 695 'Lead_Orientation
', obj.Lead_Orientation, ... 701 meta_data = metaclass(obj); 703 for idx = 1:length(meta_data.PropertyList) 704 prop_name = meta_data.PropertyList(idx).Name; 705 ret.Write( prop_name, obj.(prop_name)); 711 %> @brief Resets all elements in the object. 712 function Reset( obj ) 715 meta_data = metaclass(obj); 717 for idx = 1:length(meta_data.PropertyList) 718 prop_name = meta_data.PropertyList(idx).Name; 719 if strcmp(prop_name, 'Opt
') || strcmp( prop_name, 'param
') || strcmp(prop_name, 'varargin
') 722 obj.Clear( prop_name ); 734 %> @brief Sets the value of an attribute in the interface. 735 %> @param MemberName The name of the attribute to be set. 736 %> @param input The value to be set. 737 function Write(obj, MemberName, input) 739 if strcmp(MemberName, 'param
') 743 elseif strcmpi(MemberName, 'params
') 745 if isempty(obj.Hanyadik_Lead) 746 obj.param.scatter = input; 748 obj.param.Leads{obj.Hanyadik_Lead} = input; 752 obj.(MemberName) = input; 754 error(['EQuUs:
', class(obj), ':Read
'], ['No
property to write with name:
', MemberName]); 760 %> @brief Query for the value of an attribute in the interface. 761 %> @param MemberName The name of the attribute to be set. 762 %> @return Returns with the value of the attribute. 763 function ret = Read(obj, MemberName) 766 ret = obj.(MemberName); 768 error(['EQuUs:
', class(obj), ':Read
'], ['No
property to read with name:
', MemberName]); 773 %> @brief Clears the value of an attribute in the interface. 774 %> @param MemberName The name of the attribute to be cleared. 775 function Clear(obj, MemberName) 778 obj.(MemberName) = []; 780 error(['EQuUs:
', class(obj), ':Clear
'], ['No
property to clear with name:
', MemberName]); 788 methods ( Access = protected ) 791 %> @brief Updates the number of sites in the cross section. 793 if isempty( obj.Hanyadik_Lead ) 794 obj.M = obj.param.scatter.shape.width; 796 obj.M = obj.param.Leads{obj.Hanyadik_Lead}.M; 802 %> @brief Initializes object properties. 803 function Initialize(obj) 804 obj.InputParsing( obj.varargin{:}); 807 obj.HamiltoniansCreated = false; 808 obj.HamiltoniansDecimated = false; 809 obj.OverlapApplied = false; 810 obj.MagneticFieldApplied = false; 811 obj.GaugeTransformationApplied = false; 817 if isempty( obj.Hanyadik_Lead ) 818 obj.params = obj.param.scatter; %Lead parameters 820 obj.params = obj.param.Leads{obj.Hanyadik_Lead}; %Lead parameters 823 obj.End_Type = obj.params.End_Type; 826 end % protected methods 829 methods (Access=protected) 833 %> @brief Parses the optional parameters for the class constructor. 834 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html): 835 %> @param 'Hanyadik_Lead
' The ID number of the current lead. Set to empty (default) for using parameters of the scatter region. 836 %> @param 'Lead_Orientation
' Orientation of the lead. Set +1 (default) is the "incoming" direction of the propagating states is defined in the +x or +y direction, and "-1" otherwise. 837 %> @param 'q
' The transverse momentum. Set to empty (default) for computations without transverse momentums. 838 function InputParsing( obj, varargin ) 840 p.addParameter('Hanyadik_Lead
', []); 841 p.addParameter('Lead_Orientation
', 1); 842 p.addParameter('q
', []); 844 % keeping unmatched attributes that possibly comes from the derived classes 845 p.KeepUnmatched = true; 847 p.parse(varargin{:}); 849 obj.Lead_Orientation = p.Results.Lead_Orientation; 850 obj.Hanyadik_Lead = p.Results.Hanyadik_Lead; 858 end % private methods A class containing methodes for displaying several standard messages.
Structure Opt contains the basic computational parameters used in EQuUs.
Class to create and store Hamiltonian of the translational invariant leads.
Class to create the Hamiltonian of one unit cell in a translational invariant lead made of square lat...
function Transport(Energy, B)
creating the Ribbon class representing the twoterminal setup
function Hamiltonians(varargin)
Function to create the custom Hamiltonians for the 1D chain.
Structure containing the physical parameters describing a given lead.
Structure param contains data structures describing the physical parameters of the scattering center ...
Structure sites contains data to identify the individual sites in a matrix.
A class responsible for the Peierls and gauge transformations.
Structure containg datat on the calculated eigenstate in a translational invariant lead...
A class to create the Hamiltonian of one unit cell in a translational invariant lead made of hexagona...
A class to create and store Hamiltonian of the scattering region.