2 % Copyright (C) 2009-2015 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 A
class to evaluate the Dyson equation and to calculate the scattering matrix
for equilibrium processes.
22 %> @brief A
class to evaluate the Dyson equation and to calculate the scattering matrix for equilibrium processes.
27 properties ( Access =
protected )
28 %> The energy to be used in the calculations.
30 %> An instance of the structure
param 32 %> The calculated Greens
function.
36 %> The calculated scattering matrix.
38 %> Numbers of the open channels in the leads.
42 %> The matrix of the calculated conductance.
44 %> A
function handle of the Dyson equation.
48 %> An instance of
class #
Peierls.
50 %> An array of classes #
Lead (or its subclasses)
52 %> list of optional parameters (see http:
61 %% Contructor of the
class 62 %> @brief Constructor of the
class.
66 %> @
return An instance of the
class 72 obj.varargin = varargin;
80 %> @brief Calculates the effective (decimated) Hamiltonian of the scattering region. (Obsolete)
81 function ScatterCalc( obj )
83 obj.display(
'Creating Hamiltonian for Scattering region')
84 if isempty(obj.CreateH)
87 obj.CreateH.CreateScatterH();
89 if obj.Opt.magnetic_field == 1
90 obj.
display('Applying magnetic field in scattering region')
99 if ~isempty(obj.param.scatter.Overlap_in_Scatter) && obj.param.scatter.Overlap_in_Scatter && ~(isempty(obj.Opt.type_of_scanning) ) && ~obj.Opt.Just_Create_Hamiltonians
100 obj.
display('Applying overlap integrals in scattering region')
101 obj.ApplyOverlapInScatter( obj.CreateH, obj.E )
105 if obj.Opt.
Decimation >= 1 && ~obj.Opt.Decimate_only_Dyson && ~obj.Opt.Just_Create_Hamiltonians
107 obj.
display('Applying decimation procedure in scattering region')
109 Decimation_Procedure.DecimationFunc(obj.E, obj.CreateH, 'Hscatter', 'kulso_szabfokok');
110 obj.
display(['A decimalasnal eltelt ido:', num2str(cputime()-ido)]);
112 obj.
display('No decimation calculation is performed in the scattering region')
118 %> @brief Invokes the function
#SurfaceGreenFunctionCalculator for each lead. 119 %> @param varargin Cell array of optional parameters (https:
120 %> @param
'shiftLeads' A real number. If given, the on-site potentials of the
sites in the lead are shifted by
this value.
122 %> @param
'coordinates_shift' An integer. If given, the coordinates of the
sites in the lead are shifted by coordinates_shift*lattice vector.
123 %> @param
'gauge_field' Function
handle S = f( x,y) of the
gauge transformation on the
Hamiltonians of the leads. (S is a N x 1 vector, where N is the number of the points given by the x and y coordinates.)
124 %> @param
'createCore' Set 1
for creating the
class instances #
Lead but omitting any further calculations, or false (
default) to
continue with the calculations.
125 %> @param
'SelfEnergy' Logical value. Set
true to calculate the
self energies of the leads. (
default is
false)
126 %> @param
'SurfaceGreensFunction' Logical value. Set
true to calculate the surface Green functions of the leads. (
default is
true)
127 %> @param
'leads' Array of lead identification numbers
for which the calculations should be performed. (Default is a list of all leads)
130 %> @param
'Just_Create_Hamiltonians' Logical value. Set
true to create the
Hamiltonians of the leads, but stop further calculations. Default value is
false.
131 %> @param
'q' The tranverse momentum
for transverse computations.
132 %> @
return Returns with an array of the created #
Lead instances.
133 function lead_ret = LeadCalc( obj, varargin )
136 p.addParameter(
'shiftLeads', zeros(length(obj.param.Leads),1));
137 p.addParameter(
'transversepotential', []);
138 p.addParameter(
'coordinates_shift', zeros(length(obj.param.Leads),1) );
139 p.addParameter(
'gauge_field', [] );
140 p.addParameter(
'createCore', 0);
141 p.addParameter(
'SelfEnergy',
false); %
set true to calculate the
self energy of the semi-infinite lead
142 p.addParameter(
'SurfaceGreensFunction',
true );%
set true to calculate the surface Greens
function of the semi-infinite lead
143 p.addParameter(
'leads', 1:length(obj.param.Leads) );% list of leads to be calculated
144 p.addParameter(
'leadmodel', []); %
function handle for an individual physical model
for the contacts
145 p.addParameter(
'CustomHamiltonian', []);
146 p.addParameter(
'Just_Create_Hamiltonians', 0);
147 p.addParameter(
'q', [] ) %transverse momentum
148 p.parse(varargin{:});
150 shiftLeads = p.Results.shiftLeads;
151 transversepotential = p.Results.transversepotential;
152 coordinates_shift = p.Results.coordinates_shift; % shifts coordinates by
"coordinates_shift" times of a lattice vector.
154 createCore = p.Results.createCore;
155 SelfEnergy = p.Results.SelfEnergy;
156 SurfaceGreensFunction = p.Results.SurfaceGreensFunction;
157 leads = p.Results.leads;
158 Just_Create_Hamiltonians = p.Results.Just_Create_Hamiltonians;
160 CustomHamiltonian = p.Results.CustomHamiltonian;
161 leadmodel = p.Results.leadmodel;
163 if isempty(obj.Leads)
164 obj.Leads = cell(length(obj.param.Leads),1);
167 obj.display(
'Creating surface_Greens_function interfaces:')
170 while idx <= length(leads)
171 obj.display(
'------------------------------')
172 obj.display([
'Lead: ',num2str(leads(idx))])
173 obj.Leads{leads(idx)} = obj.SurfaceGreenFunctionCalculator(leads(idx),
'shiftLead', shiftLeads(leads(idx)),
'transversepotential', transversepotential,
'coordinates_shift', coordinates_shift(leads(idx)),...
174 'gauge_field',
gauge_field,
'createCore', createCore, ...
175 'SelfEnergy', SelfEnergy,
'SurfaceGreensFunction', SurfaceGreensFunction,
'q', q,
'leadmodel', leadmodel,
'CustomHamiltonian', CustomHamiltonian,
'Just_Create_Hamiltonians', Just_Create_Hamiltonians);
180 lead_ret = obj.Leads;
185 %> @brief Invokes the
function handle of the Dyson equation, and stores the calculated values in attributes #G and #
junction_sites.
186 %> @param varargin Cell array of optional parameters (https:
187 %> @param
'CustomDyson' A
function handle [G, Ginverz, #
junction_sites] = CustomDyson() of the Dyson equation to be evaluated. If not given, the
function handle given by the present
class is used instead.
188 %> @
return [1] Returns with the calculated Green
operator.
190 function [Gret, junction_sites_ret] = DysonEq( obj, varargin )
192 obj.display([
'EQuUs:',
class(obj),
'DysonEq: Evaluating Dyson Equation'])
195 p.addParameter(
'CustomDyson', []);
196 p.parse(varargin{:});
198 if ~isempty( p.Results.CustomDyson )
199 obj.Dysonfunc = p.Results.CustomDyson;
202 if ~isempty( obj.Dysonfunc )
203 if nargout( obj.Dysonfunc ) == 1
204 obj.G = obj.Dysonfunc();
206 [obj.G, ~, obj.junction_sites] = obj.Dysonfunc();
212 junction_sites_ret = obj.junction_sites;
218 %> @brief Conductance calculated by Eq (19) in PRB 73 085414
219 %> @
return Returns with the calculated conductance.
220 function C = Conductance2( obj )
223 Gamma = cell(size(obj.Leads));
224 for idx = 1:length(obj.Leads)
225 Gamma{idx} = obj.Leads{idx}.Gamma();
228 Gamma_left = [Gamma{1}, zeros( size(Gamma{1},1), size(Gamma{2},2));
229 zeros( size(Gamma{2},1), size(Gamma{1},2)+size(Gamma{2},2)) ];
231 Gamma_right = [zeros( size(Gamma{1},1), size(Gamma{1},2)+size(Gamma{2},2));
232 zeros(size(Gamma{2},1), size(Gamma{1},2)), Gamma{2}];
234 for idx = 1:length(obj.Leads)
238 % Eq (19) PRB 73 085414
239 C = trace( Gamma_left*obj.G
'*Gamma_right*obj.G ) ; 242 %*********************************** 243 % Sigma = cell(size(obj.Leads)); 244 % for idx = 1:length(obj.Leads) 245 % Sigma{idx} = obj.Leads{idx}.Read('Sigma
'); 248 % Sigma_left = [Sigma{1}, zeros( size(Sigma{1},1), size(Sigma{2},2)); 249 % zeros( size(Sigma{2},1), size(Sigma{1},2)+size(Sigma{2},2)) ]; 251 % Sigma_right = [zeros( size(Sigma{1},1), size(Sigma{1},2)+size(Sigma{2},2)); 252 % zeros(size(Sigma{2},1), size(Sigma{1},2)), Sigma{2}]; 254 % for idx = 1:length(obj.Leads) 258 % Eq (22) Eur. Phys. J. B 53, 537-549 (2006) 259 % current2 = 2*real( trace(-1i*obj.G*Gamma_right*obj.G'*Sigma_left
') ); 262 %********************************* 268 %> @brief Calculates the conductance matrix from the scattering matrix 269 %> @return Returns with the calculated conductance matrix. 270 function C = Conduktance( obj ) 271 obj.display('calculating conductance matrix
') 273 obj.nyitott_csatornak = obj.nyitott_csatornak; 274 LeadsNumber = length(obj.open_channels); 275 C = zeros(LeadsNumber, LeadsNumber); 278 for rowidx = 1:LeadsNumber 279 rows_tmp = max(rows)+ (1:obj.open_channels(rowidx).num); 285 for colidx = 1:LeadsNumber 286 cols_tmp = max(cols)+ (1:obj.open_channels(colidx).num); 291 s = obj.Smatrix( rows, cols); 296 C(rowidx,colidx) = -obj.open_channels(rowidx).num + trace(s'*s);
298 C(rowidx,colidx) = trace(s
'*s); 301 if obj.open_channels(rowidx).isSuperconducting || obj.open_channels(colidx).isSuperconducting 304 incoming_electrons = obj.open_channels(colidx).BdG_u_incoming( obj.open_channels(colidx).open_channels_p ); 305 outgoing_holes = ~obj.open_channels(rowidx).BdG_u_outgoing( obj.open_channels(colidx).open_channels_p ); 306 outgoing_electrons = obj.open_channels(rowidx).BdG_u_outgoing( obj.open_channels(colidx).open_channels_m ); 309 r_normal = s( outgoing_electrons, incoming_electrons ); 310 r_andreev = s( outgoing_holes, incoming_electrons ); 311 M = length(find( obj.open_channels(rowidx).open_channels_p & obj.open_channels(rowidx).BdG_u_incoming)); 312 C(rowidx,colidx) = -M + trace(r_normal'*r_normal) - trace(r_andreev
'*r_andreev); 314 t_normal = s( outgoing_electrons, incoming_electrons ); 315 t_andreev = s( outgoing_holes, incoming_electrons ); 316 C(rowidx,colidx) = trace(t_normal'*t_normal) - trace(t_andreev
'*t_andreev); 324 obj.Conductance_Matrix = C; 330 %> @brief Calculates the scattering matrix 331 %> @return [1] The scattering matrix 332 %> @return [2] An array containing the number of openchannels in the leads. Obsolete, use attribute #open_channels istead. 333 function [S, Neff] = SmatrixCalc( obj ) 334 obj.display('calculating S-matrix
') 335 Neff = obj.Get_Neff(); 336 obj.Determine_Open_Channels() 338 NormamtxAll = zeros(sum(Neff), sum(Neff)); 339 Modusmtx_pAll = zeros(sum(Neff), sum(Neff)); 340 d_Modusmtx_mAll = zeros(sum(Neff), sum(Neff)); 341 Gyokvcsop_pAll = zeros(sum(Neff), 1); 342 Gyokvcsop_mAll = zeros(sum(Neff), 1); 343 evanescent_indexes_pAll = false(sum(Neff), 1); 344 evanescent_indexes_mAll = false(sum(Neff), 1); 346 for idx = 1:length(Neff) 347 NormamtxAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)), sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = obj.Leads{idx}.Read('Normamtx
'); 348 Modusmtx_pAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)), sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = obj.Leads{idx}.Read('modusmtx_p
'); 349 d_Modusmtx_mAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)), sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = obj.Leads{idx}.Read('d_modusmtx_m
'); 350 Gyokvcsop_pAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = sqrt(abs(obj.Leads{idx}.Read('vcsop_p
'))); 351 Gyokvcsop_mAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = sqrt(abs(obj.Leads{idx}.Read('vcsop_m
'))); 352 evanescent_indexes_pAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = ~(obj.open_channels(idx).open_channels_p); 353 evanescent_indexes_mAll( sum(Neff(1:idx-1))+1:sum(Neff(1:idx)) ) = ~(obj.open_channels(idx).open_channels_m); 356 GN0 = obj.G(1:sum(Neff),1:sum(Neff))*NormamtxAll; 357 S = (d_Modusmtx_mAll*(GN0-eye(sum(Neff))) )*Modusmtx_pAll; 359 % sorting out evanescent modes 360 S(:,evanescent_indexes_pAll) = []; 361 S(evanescent_indexes_mAll,:) = []; 362 Gyokvcsop_mAll = Gyokvcsop_mAll(~evanescent_indexes_mAll); 363 Gyokvcsop_pAll = Gyokvcsop_pAll(~evanescent_indexes_pAll); 365 S = diag( Gyokvcsop_mAll )*S*diag(Gyokvcsop_pAll.^(-1)); 368 % obsolete code lines 369 obj.nyitott_csatornak = zeros(size(Neff)); 370 for idx = 1:length(Neff) 371 obj.nyitott_csatornak(idx) = obj.open_channels(idx).num; 373 Neff = obj.nyitott_csatornak; 379 %% Determine_Open_Channels 380 %> @brief Determines the open channels in the leads. The data are storen within the attribute #open_channels. 381 function Determine_Open_Channels( obj ) 384 obj.open_channels = repmat(obj.open_channels, size(M)); 386 for idx = 1:length( M ) 387 obj.Leads{idx}.Determine_Open_Channels(); 388 obj.open_channels(idx) = obj.Leads{idx}.Read('open_channels'); 393 %% SurfaceGreenFunctionCalculator 394 %> @brief Calculates the surface Green's
function or the
self energy of a
Lead.
396 %> @param varargin Cell array of optional parameters (https:
397 %> @param
'shiftLeads' A real number. If given, the on-site potentials of the
sites in the lead are shifted by
this value.
399 %> @param
'Lead' An instance of
class #
Lead (or its subclass) that is intended to be used in the calculations.
400 %> @param
'coordinates_shift' An integer. If given, the coordinates of the
sites in the lead are shifted by coordinates_shift*lattice vector.
401 %> @param
'gauge_field' Function
handle S = f( x,y) of the
gauge transformation on the
Hamiltonians of the leads. (S is a N x 1 vector, where N is the number of the points given by the x and y coordinates.)
402 %> @param
'createCore' Set 1
for creating the
class instances #
Lead but omitting any further calculations, or false (
default) to
continue with the calculations.
403 %> @param
'SelfEnergy' Logical value. Set
true to calculate the
self energies of the leads. (
default is
false)
404 %> @param
'SurfaceGreensFunction' Logical value. Set
true to calculate the surface Green functions of the leads. (
default is
true)
405 %> @param
'leads' Array of lead identification numbers
for which the calculations should be performed. (Default is a list of all leads)
408 %> @param
'Just_Create_Hamiltonians' Logical value. Set
true to create the
Hamiltonians of the leads, but stop further calculations. Default value is
false.
409 %> @param
'q' The tranverse momentum
for transverse computations.
410 %> @
return Returns with the created #
Lead instance.
411 function Lead_ret = SurfaceGreenFunctionCalculator( obj, idx, varargin)
414 p.addParameter(
'createCore', 0);
415 p.addParameter(
'Just_Create_Hamiltonians', 0);
416 p.addParameter(
'shiftLead', 0);
417 p.addParameter(
'coordinates_shift', 0);
418 p.addParameter(
'transversepotential',[]);
419 p.addParameter(
'Lead',[]);
420 p.addParameter(
'gauge_field', [] );%
gauge field
for performing
gauge transformation
421 p.addParameter(
'SelfEnergy',
false); %
set true to calculate the
self energy of the semi-infinite lead
422 p.addParameter(
'SurfaceGreensFunction',
true );%
set true to calculate the surface Greens
function of the semi-infinite lead
423 p.addParameter(
'leadmodel', []); %
function handle for an individual physical model
for the contacts
424 p.addParameter(
'CustomHamiltonian', []);
425 p.addParameter(
'q', [] ) %transverse momentum
426 p.parse(varargin{:});
427 createCore = p.Results.createCore;
428 Just_Create_Hamiltonians = p.Results.Just_Create_Hamiltonians;
429 shiftLead = p.Results.shiftLead;
430 coordinates_shift = p.Results.coordinates_shift;
431 transversepotential = p.Results.transversepotential;
432 Lead_ret = p.Results.Lead;
434 SelfEnergy = p.Results.SelfEnergy;
435 SurfaceGreensFunction = p.Results.SurfaceGreensFunction;
437 CustomHamiltonian = p.Results.CustomHamiltonian;
438 leadmodel = p.Results.leadmodel;
441 if ~isempty( Lead_ret )
442 %
do not create a
new class instance if Lead_ret was given
445 elseif ~isempty(idx) && ~isempty(obj.param.Leads{idx}.Lead_Orientation)
446 Lead_ret =
Lead(obj.Opt, obj.param,...
447 'Hanyadik_Lead', idx,....
448 'Lead_Orientation', obj.param.Leads{idx}.Lead_Orientation, ...
451 Lead_ret =
Lead(obj.Opt, obj.param,...
452 'Hanyadik_Lead', idx, ...
461 if obj.Opt.Simple_Green_Function
462 Lead_ret.SurfaceGreen_simple(obj.E);
463 elseif ~isempty(leadmodel)
464 Lead_ret = leadmodel( idx, obj.E,
'createCore', createCore, ...
465 'Just_Create_Hamiltonians', Just_Create_Hamiltonians, ...
466 'shiftLead', shiftLead, ...
467 'coordinates_shift', coordinates_shift, ...
468 'transversepotential', transversepotential, ...
469 'Lead', Lead_ret, ...
471 'SelfEnergy', SelfEnergy, ...
472 'SurfaceGreensFunction', SurfaceGreensFunction, ...
476 % creating Hailtonians
477 obj.display([
'EQuUs:',
class(obj),
':SurfaceGreenFunctionCalculator: Creating Hamiltonians for lead'])
478 Lead_ret.CreateHamiltonians(
'CustomHamiltonian', CustomHamiltonian );
479 Lead_ret.ShiftCoordinates( coordinates_shift );
480 if obj.Opt.magnetic_field == 1
481 obj.display([
'EQuUs:',
class(obj),
':SurfaceGreenFunctionCalculator: Applying magnetic field in lead'])
483 % In superconducting lead one must not include nonzero magnetic
486 % traslational invariant. In transverse computations a gauge
487 % tranformation is performed on each lead-sample
interface to
488 % cancel the vector potential
489 if ~Lead_ret.isSuperconducting() %&& isempty( Lead_tmp.Read(
'q') )
490 obj.PeierlsTransform.PeierlsTransformLeads(Lead_ret);
491 elseif ~isempty(
gauge_field ) %&& isempty( Lead_tmp.Read(
'q') )
492 obj.PeierlsTransform.gaugeTransformationOnLead( Lead_ret,
gauge_field );
496 if abs(shiftLead) > 1e-6
497 Lead_ret.ShiftLead( shiftLead );
500 if ~isempty(transversepotential) && isempty(q) %In transverse computations no transverse potential can be applied
501 coordinates = Lead_ret.Read(
'coordinates');
502 if nargin( transversepotential ) == 1
503 potential2apply = transversepotential( coordinates );
504 elseif nargin( transversepotential ) == 2
505 potential2apply = transversepotential( Lead_ret, obj.E );
507 error(
'EQuUs:Transport_Interface:SurfaceGreenFunctionCalculator',
'To many input arguments in function handle PotInScatter');
510 if isprop(coordinates,
'BdG_u')
511 if size( potential2apply, 1) == 1 || size( potential2apply, 2) == 1
512 potential2apply(~coordinates.BdG_u) = -potential2apply(~coordinates.BdG_u);
514 potential2apply(~coordinates.BdG_u, ~coordinates.BdG_u) = -conj(potential2apply(~coordinates.BdG_u, ~coordinates.BdG_u));
517 Lead_ret.AddPotential( potential2apply );
520 if Just_Create_Hamiltonians
524 % Solve the eigen problem
525 obj.display([
'EQuUs:',
class(obj),
':SurfaceGreenFunctionCalculator: Eigenvalues for lead'])
526 Lead_ret.TrukkosSajatertekek(obj.E);
529 obj.display([
'EQuUs:',
class(obj),
':SurfaceGreenFunctionCalculator: Group velocities for lead'])
530 Lead_ret.Group_Velocity();
533 % retarded surface Green
operator 534 if SurfaceGreensFunction
535 Lead_ret.SurfaceGreenFunction();
538 % retarded SelfEnergy
540 Lead_ret.SelfEnergy();
549 %> @brief Sets the energy to be used in the calculations and resets the calculated attributes.
550 %> @param Energy The energy value
551 function setEnergy( obj, Energy )
555 if ~isempty( obj.CreateH )
559 % reset the lead attributes
560 if ~isempty( obj.Leads )
561 for idx = 1:length( obj.Leads )
562 obj.Leads{idx}.Reset();
570 %> @brief Query
for the value of an attribute in the
class.
571 %> @param MemberName The name of the attribute.
572 %> @
return Returns with the value of the attribute.
573 function ret = Read(obj, MemberName)
574 if strcmp(MemberName,
'M')
579 ret = obj.(MemberName);
581 error([
'EQuUs:',
class(obj),
':Read'], [
'No property to read with name: ', MemberName]);
588 %> @brief Sets the value of an attribute in the
class.
589 %> @param MemberName The name of the attribute to be
set.
590 %> @param input The value to be
set.
591 function Write(obj, MemberName, input)
593 if strcmp(MemberName,
'param')
597 elseif strcmp(MemberName,
'E')
598 obj.setEnergy( input );
601 obj.(MemberName) = input;
603 error([
'EQuUs:',
class(obj),
':Read'], [
'No property to write with name: ', MemberName]);
610 %> @brief Clears the value of an attribute in the
class.
611 %> @param MemberName The name of the attribute to be cleared.
612 function Clear(obj, MemberName)
614 obj.(MemberName) = [];
616 error([
'EQuUs:',
class(obj),
':Clear'], [
'No property to clear with name: ', MemberName]);
622 %> @brief Resets all elements in the
class.
625 if strcmpi(
class(obj),
'Transport_Interface' )
626 meta_data = metaclass(obj);
628 for idx = 1:length(meta_data.PropertyList)
629 prop_name = meta_data.PropertyList(idx).Name;
630 if strcmp(prop_name,
'Opt') || strcmp( prop_name,
'param') || strcmp(prop_name,
'varargin') || strcmp(prop_name,
'E')
633 obj.Clear( prop_name );
644 %> @brief Adds/replaces a lead to/in the system.
645 %> @param Lead_tmp An instance of
class #
Lead 647 function replaceLead( obj, Lead_tmp, idx )
649 if ~strcmpi(
class(obj.junction),
'Lead' )
650 supClasses = superclasses(Lead_tmp);
651 if sum( strcmp( supClasses,
'Lead' ) ) == 0
652 error([
'EQuUs:',
class(obj),
':replaceLead'],
'Input Lead_tmp is not valid.');
656 obj.Leads{ idx } = Lead_tmp;
657 obj.Opt.NofLeads = length(obj.Leads);
662 function LeadSave( obj )
663 for idx = 1:length(obj.Leads)
664 obj.Leads{idx}.SaveLeads();
669 %> @brief Creates a clone of the present
class.
670 %> @
return Returns with the cloned
object.
673 if isempty( obj.CreateH )
676 CreateH_loc = obj.CreateH.CreateClone();
679 if isempty( obj.PeierlsTransform )
680 PeierlsTransform_loc = [];
682 PeierlsTransform_loc = obj.PeierlsTransform.CreateClone();
686 'CreateH', CreateH_loc, ...
687 'PeierlsTransform', PeierlsTransform_loc );
689 for idx = 1:length(obj.Leads)
690 ret.replaceLead( obj.Leads{idx}.CreateClone(), idx );
697 %> @
return Returns with an array of the lead widths
698 function M = GetM( obj )
699 M = zeros(length(obj.Leads),1);
700 for idx = 1:length(obj.Leads)
701 M(idx) = obj.Leads{idx}.Read(
'M' );
710 %> @brief Gets the size of the effective
Hamiltonians of the leads
711 %> @
return Returns with an array of integers;
713 Neffs = zeros(length(obj.Leads),1);
714 for idx = 1:length(obj.Leads)
715 Neffs(idx) = obj.Leads{idx}.Get_Neff();
723 methods ( Access =
protected )
726 %> @brief Initializes
object attributes.
730 obj.junction_sites = [];
732 obj.nyitott_csatornak = [];
733 obj.open_channels = [];
734 obj.Conductance_Matrix = [];
738 obj.PeierlsTransform = [];
741 if strcmpi(
class(obj),
'Transport_Interface')
742 obj.InputParsing( obj.varargin{:});
748 %> @brief Parses the optional parameters
for the
class constructor.
749 %> @param varargin Cell array of optional parameters (https:
750 %> @param
'CreateH' An instance of
class #
CreateHamiltonians storing the Hamiltonian of the scattering center.
751 %> @param
'PeierlsTransform' An instance of
class #
Peierls for the
Peierls transformations.
752 function InputParsing( obj, varargin )
755 p.addParameter(
'CreateH', []);
756 p.addParameter(
'PeierlsTransform', []);
758 p.parse(varargin{:});
761 obj.CreateH = p.Results.CreateH;
765 end % methods
protected function Initialize()
Eq (20) in PRB 78, 035407 (2008)
function Messages(Opt)
Constructor of the class.
function InputParsing(varargin)
Parses the optional parameters for the class constructor.
Property Hanyadik_Lead
The id number of the current lead.
A class containing methodes for displaying several standard messages.
Structure Opt contains the basic computational parameters used in EQuUs.
Property gauge
String containing the name of the built-in gauge ('LandauX', 'LandauY').
Structure open_channels describes the open channel in the individual leads.
Class to create and store Hamiltonian of the translational invariant leads.
function display(message, nosilent)
Displays output messages on the screen.
function Transport(Energy, B)
creating the Ribbon class representing the twoterminal setup
Property M
The number of the sites in the cross section.
function PeierlsTransform(CreateH)
Algorithm to perform Peierls transformation in Hamiltonian stored in the CreateHamiltonians object...
A class to evaluate the Dyson equation and to calculate the scattering matrix for equilibrium process...
function Hamiltonians(varargin)
Function to create the custom Hamiltonians for the 1D chain.
A class providing function handle to reduce the number of sites in the Hamiltonian via decimation pro...
function InputParsing(varargin)
Parses the optional parameters for the class constructor.
A class to import custom Hamiltonians provided by other codes or created manually ...
function Reset()
Resets all elements in the class.
A class to calculate the Green functions and self energies of a translational invariant lead The nota...
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.
function Get_Neff()
Gets the effective number of sites after the elimination of the singular values.
function SurfaceGreenFunctionCalculator(idx, varargin)
Calculates the surface Green's function or the self energy of a Lead.
A class responsible for the Peierls and gauge transformations.
function gauge_field(x, y, eta_B, Aconst, height, lattice_constant)
Scalar gauge field connecting Landaux and Landauy gauges.
function CreateClone()
Creates a clone of the present object.
Structure junction_sites contains data to identify the individual sites of the leads, the interface regions and the scattering region in a matrix.
A class to create and store Hamiltonian of the scattering region.