Eötvös Quantum Utilities  v4.8.141
Providing the Horsepowers in the Quantum Realm
Ribbon.m
Go to the documentation of this file.
1 %% Eotvos Quantum Transport Utilities - Ribbon
2 % Copyright (C) 2009-2015 Peter Rakyta, Ph.D.
3 %
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.
8 %
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.
13 %
14 % You should have received a copy of the GNU General Public License
15 % along with this program. If not, see http://www.gnu.org/licenses/.
16 %
17 %> @addtogroup utilities Utilities
18 %> @{
19 %> @file Ribbon.m
20 %> @brief A class for calculations on a ribbon of finite width for equilibrium calculations mostly in the zero temperature limit.
21 %> @image html Ribbon_structure.jpg
22 %> @image latex Ribbon_structure.jpg
23 %> @}
24 %> @brief A class for calculations on a ribbon of finite width for equilibrium calculations mostly in the zero temperature limit.
25 %> <tr class="heading"><td colspan="2"><h2 class="groupheader"><a name="avail"></a>Structure described by the class</h2></td></tr>
26 %> @image html Ribbon_structure.jpg
27 %> @image latex Ribbon_structure.jpg
28 %> The drawing represents a two-terminal structure made of two leads, of a scattering region and of two interface regions between the leads and scattering center.
29 %> Each rectangle describes a unit cell including singular and non-singular sites.
30 %> The scattering center is also described by a set of identical unit cells, but arbitrary potential can be used.
31 %> Arrows indicate the hopping direction stored by the attributes in the corresponding classes (see attributes #CreateLeadHamiltonians.H1 and #InterfaceRegion.Hcoupling for details).
32 %> The orientation of the lead is +1 if the lead is terminated by the interface region in the positive direction, and -1 if the lead is terminated by the interface region in the negative direction.
33 %> (see attribute #CreateLeadHamiltonians.Lead_Orientation for details)
34 classdef Ribbon < NTerminal
35 
36 
37  properties ( Access = public )
38  %> An instance of class #Lead (or its subclass) describing the unit cell of the scattering region
39  Scatter_UC
40  %> A function handle pot = f( #coordinates ) or pot=f( #CreateLeadHamiltonians, Energy) of the transverse potential applied in the lead. (Instead of #CreateLeadHamiltonians can be used its derived class)
41  transversepotential
42  %> width of the scattering region (number of the nonsingular atomic sites in the cross section)
43  width
44  %> height (length) of the scattering region (number of unit cells)
45  height
46  %> the shift of the coordinates of the sites (two component vector)
47  shift
48  end
49 
50 
51 methods ( Access = public )
52 %% Contructor of the class
53 %> @brief Constructor of the class.
54 %> @param Opt An instance of the structure #Opt.
55 %> @param varargin Cell array of optional parameters. For details see #InputParsing.
56 %> @return An instance of the class
57  function obj = Ribbon( varargin )
58  obj = obj@NTerminal();
59 
60 
61  obj.param = [];
62  obj.Scatter_UC = [];
63  obj.transversepotential = [];
64  obj.width = [];
65  obj.height = [];
66  obj.shift = 0;
67 
68 
69  if strcmpi( class(obj), 'Ribbon')
70  % processing the optional parameters
71  obj.InputParsing(varargin{:})
72 
73 
74  obj.display(['EQuUs:Utils:', class(obj), ':Ribbon: Creating a Ribbon object'])
75 
76  % create the shape of the scattering region
77  obj.createShape();
78 
79  % set the Fermi level
80  obj.setFermiEnergy();
81 
82  % exporting calculation parameters into an XML format
83  createOutput( obj.filenameOut, obj.Opt, obj.param );
84 
85  % create class intances and initializing class attributes
86  obj.CreateHandles();
87 
88  % create Hamiltonians and coordinates of the unit cells
89  obj.CreateRibbon('justHamiltonians', true);
90  end
91 
92 
93 
94  end
95 
96 %% Transport
97 %> @brief Calculates the transport through the two terminal setup on two dimensional lattices. Use for development pupose only.
98 %> @param Energy The energy value.
99 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
100 %> @param 'constant_channels' Logical value. Set true (default) to keep constant the number of the open channels in the leads for each energy value, or false otherwise.
101 %> @param 'gfininvfromHamiltonian' Logical value. Set true calculate the surface Greens function of the scattering region from the Hamiltonaian of the scattering region, or false (default) to calculate it by the fast way (see Phys. Rev. B 90, 125428 (2014) for details).
102 %> @param 'decimateDyson' Logical value. Set true (default) to decimate the sites of the scattering region in the Dyson equation.
103 %> @param 'PotInScatter' Obsolete parameter. Use 'scatterPotential' instead.
104 %> @param 'scatterPotential' A function handle pot=f( #coordinates ) or pot=f( #CreateHamiltonians, Energy) for the potential to be applied in the Hamiltonian (used when FiniteGreensFunctionFromHamiltonian=true).
105 %> @param 'selfEnergy' Logical value. Set true to use the self energies of the leads in the Dyson equation, or false (default) to use the surface Green function instead.
106 %> @param 'Smatrix' Set true (default) to calculate the conductance by using the scattering matrix via #Transport_Interface.Conduktance, or false to use function #Transport_Interface.Conductance2. (The latter one also works with complex energies.)
107 %> @return [1] Conductivity The calculated conductivity.
108 %> @return [2] aspect_ratio The aspect ratio W/L of the junction.
109 %> @return [3] Conductance The conductance tensor
110 %> @return [4] ny Array of the open channel in the leads.
111 %> @return [5] DeltaC Error of the unitarity.
112 %> @return [6] S The scattering matrix.
113  function [Conductivity,aspect_ratio,Conductance,ny,DeltaC,S] = Transport(obj, Energy, varargin)
114 
115  p = inputParser;
116  p.addParameter('constant_channels', false);
117  p.addParameter('gfininvfromHamiltonian', false);
118  p.addParameter('decimateDyson', true);
119  p.addParameter('PotInScatter', []) %OBSOLETE use scatterPotential instead
120  p.addParameter('scatterPotential', []) %NEW overrides optional argument 'PotInScatter'
121  p.addParameter('selfEnergy', false)
122  p.addParameter('Smatrix', true) %logcal value to use the S-matrix method to calculate the conductance or not.
123  p.parse(varargin{:});
124  constant_channels = p.Results.constant_channels;
125  gfininvfromHamiltonian = p.Results.gfininvfromHamiltonian;
126 
127  scatterPotential = p.Results.PotInScatter;
128  if ~isempty( p.Results.scatterPotential )
129  scatterPotential = p.Results.scatterPotential;
130  end
131 
132  decimateDyson = p.Results.decimateDyson;
133  selfEnergy = p.Results.selfEnergy;
134  Smatrix = p.Results.Smatrix;
135 
136  obj.CalcSpectralFunction(Energy, 'constant_channels', constant_channels, 'gfininvfromHamiltonian', gfininvfromHamiltonian, ...
137  'decimateDyson', decimateDyson, 'scatterPotential', scatterPotential, 'SelfEnergy', selfEnergy);
138 
139  if strcmp(obj.param.scatter.End_Type, 'A')
140  aspect_ratio = (obj.width*3/2)/((obj.height+3)*sqrt(3));
141  else
142  aspect_ratio = ((obj.width-0.5)*sqrt(3))/((obj.height+3)*3);
143  end
144 
145  if Smatrix
146  [S,ny] = obj.FL_handles.SmatrixCalc();
147 
148  norma = norm(S*S'-eye(sum(ny)));
149 
150  if norma >= 1e-3
151  obj.display( ['error of the unitarity of S-matrix: ',num2str(norma)] )
152  end
153 
154  if ny(1) ~= ny(2)
155  obj.display( ['openchannels do not match: ',num2str(ny(1)), ' ', num2str(ny(2))] )
156  end
157 
158 
159 
160  conductance = obj.FL_handles.Conduktance();
161  conductance = abs([conductance(1,:), conductance(2,:)]);
162  DeltaC = std(conductance);
163 
164  C = mean(conductance);
165  Conductivity = C/aspect_ratio;
166  Conductance = C;
167 
168  obj.display( ['aspect ratio = ', num2str(aspect_ratio), ' conductance = ', num2str(C), ' open_channels= ', num2str(ny(1)), ' conductivity = ', num2str(Conductivity)])
169 
170  else
171  Conductance = obj.FL_handles.Conductance2();
172  DeltaC = NaN;
173  Conductivity = Conductance/aspect_ratio;
174  ny = NaN;
175  S = [];
176  obj.display( ['aspect ratio = ', num2str(aspect_ratio), ' conductance = ', num2str(Conductance), ' conductivity = ', num2str(Conductivity)])
177  end
178 
179  end
180 
181 %% getCoordinates
182 %> @brief Gets the coordinates of the central region
183 %> @return [1] Coordinates of the central region.
184 %> @return [2] Coordinates of the interface region.
185 function [coordinates, coordinates_interface] = getCoordinates( obj )
186  try
187  coordinates = obj.Scatter_UC.Read('coordinates');
188  H1 = obj.Scatter_UC.Read('H1');
189  [rows, cols] = find(H1);
190  rows = unique(rows);
191  cols = unique(cols);
192  fnames = fieldnames(coordinates);
193  for idx = 1:length( fnames )
194  fname = fnames{idx};
195  if strcmp( fname, 'a' ) || strcmp( fname, 'b' ) || strcmp( fname, 'LatticeConstant' )
196  continue
197  elseif strcmp( fname, 'x' )
198  coordinates.(fname) = [coordinates.(fname)(cols); coordinates.(fname)(rows) + (obj.height-1)*coordinates.a(1)];
199  elseif strcmp( fname, 'y' )
200  coordinates.(fname) = [coordinates.(fname)(cols); coordinates.(fname)(rows) + (obj.height-1)*coordinates.a(2)];
201  else
202  if ~isempty(coordinates.(fname))
203  coordinates.(fname) = [coordinates.(fname)(cols); coordinates.(fname)(rows)];
204  end
205  end
206  end
207 
208  catch errCause
209  err = MException('EQuUs:Utils:Ribbon:getCoordinatesFromRibbon', 'Error occured when retriving the coordinates of the ribbon.');
210  err = addCause(err, errCause);
211  save('Error_Ribbon_getCoordinatesFromRibbon.mat');
212  throw( err );
213  end
214 
215 
216  try
217  if isempty( obj.Interface_Regions )
218  coordinates_interface = [];
219  return
220  end
221 
222  coordinates_interface = cell( length(obj.Interface_Regions), 1);
223  for idx = 1:length( obj.Interface_Regions )
224  coordinates_interface{idx} = obj.Interface_Regions{idx}.Read('coordinates');
225  end
226 
227  catch errCause
228  err = MException('EQuUs:Utils:Ribbon:getCoordinatesFromRibbon', 'Error occured when retriving the coordinates of the interface region from a Ribbon interface');
229  err = addCause(err, errCause);
230  save('Error_Ribbon_getCoordinatesFromRibbon.mat');
231  throw( err );
232  end
233 
234 end
235 
236 %% ShiftCoordinates
237 %> @brief Shifts the coordinates of the sites in the ribbon by an integer multiple of the lattice vector. The coordinates of the Leads are automatically adjusted later
238 %> @param shift An integer value.
239  function ShiftCoordinates( obj, shift )
240  obj.Scatter_UC.ShiftCoordinates( shift );
241  if ~isempty( obj.Interface_Regions )
242  for idx = 1:length(obj.Interface_Regions)
243  obj.Interface_Regions{idx}.ShiftCoordinates( shift );
244  end
245  end
246 
247  obj.shift = obj.shift + shift;
248  end
249 
250 %% CreateScatter
251 %> @brief Initializes class #CreateHamiltonians for storing and manipulate the Hamiltonian of the the scattering region. The created object is stored in attribute #CreateH.
252  function Scatter_UC = CreateScatter( obj )
253  Scatter_UC = obj.Scatter_UC.CreateClone('empty', true);
254 
255  % create Hamiltonian of one unit cell of the scattering region
256  Scatter_UC.CreateHamiltonians( 'toSave', 0);
257  Scatter_UC.ShiftCoordinates( obj.shift );
258 
259  %applying transverse potential
260  obj.ApplyTransversePotential( Scatter_UC )
261 
262 
263  % apply magnetic field in the unit cell of the scattering region
264  % can be applied if the vector potential is identical in each unit cells
265  if ~isempty( obj.PeierlsTransform_Scatter ) && isempty(obj.q) && obj.Opt.magnetic_field_trans_invariant %for finite q the vector potential must be parallel to q, and perpendicular to the unit cell vector
266  obj.display(['EQuUs:Utils:',class(obj),':CreateScatter: Applying magnetic field in the unit cell of the scattering region']);
267  obj.PeierlsTransform_Scatter.PeierlsTransformLeads( Scatter_UC );
268  end
269 
270  % Create the Hamiltonian of the scattering region
271  createH = CreateHamiltonians(obj.Opt, obj.param, 'q', obj.q);
272  createH.CreateScatterH( 'Scatter_UC', Scatter_UC );
273 
274  % apply magnetic field in the whole HAmiltonian of the scattering region
275  % can be applied for non-translational invariant vector potentials
276  if ~isempty( obj.PeierlsTransform_Scatter ) && isempty(obj.q) && ~obj.Opt.magnetic_field_trans_invariant
277  obj.display(['EQuUs:Utils:',class(obj),':CreateScatter: Applying magnetic field in the whole Hamiltonian of the scattering region']);
278  obj.PeierlsTransform_Scatter.PeierlsTransform( createH );
279  end
280 
281 
282  obj.CreateH = createH;
283 
284  H0 = Scatter_UC.Read('H0');
285  H1 = Scatter_UC.Read('H1');
286 
287  [rows, cols] = find( H1 );
288  rows = unique(rows); % cols identical to non_singular_sites
289  cols = unique(cols); % cols identical to non_singular_sites
290 
291  % identify sites that would be directly connected to the leads
292  non_singular_sites = [reshape(cols, 1, length(cols)), (obj.height-1)*size(H0,1)+reshape(rows, 1, length(rows))];
293  createH.Write('kulso_szabfokok', non_singular_sites );
294 
295 
296  end
297 
298 %% setEnergy
299 %> @brief Sets the energy for the calculations
300 %> @param Energy The value of the energy in the same units as the Hamiltonian.
301  function setEnergy( obj, Energy )
302 
303  setEnergy@NTerminal( obj, Energy );
304 
305  if ~isempty( obj.Scatter_UC ) && strcmpi(class(obj.Scatter_UC), 'Lead')
306  obj.Scatter_UC.Reset();
307  end
308 
309  if ~isempty(obj.Scatter_UC)
310  obj.CreateRibbon('justHamiltonians', true);
311  end
312 
313  end
314 
315 
316 %% CustomDysonFunc
317 %> @brief Custom Dyson function for a two terminal arrangement on a two dimensional lattice.
318 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
319 %> @param 'gfininv' The inverse of the Greens function of the scattering region. For default the inverse of the attribute #G is used.
320 %> @param 'constant_channels' Logical value. Set true (default) to keep constant the number of the open channels in the leads for each energy value, or false otherwise.
321 %> @param 'onlyGinverz' Logical value. Set true to calculate only the inverse of the total Green operator, or false (default) to calculate #G as well.
322 %> @param 'recalculateSurface' A vector of the identification numbers of the lead surfaces to be recalculated.
323 %> @param 'decimate' Logical value. Set true (default) to eliminate all inner sites in the Greens function and keep only the selected sites. Set false to omit the decimation procedure.
324 %> @param 'kulso_szabfokok' Array of sites to be kept after the decimation procedure. (Use parameter 'keep_sites' instead)
325 %> @param 'selfEnergy' Logical value. Set true to use the self energies of the leads in the Dyson equation, or false (default) to use the surface Green function instead.
326 %> @param 'keep_sites' Name of sites to be kept in the resulted Green function (Possible values are: 'scatter', 'interface', 'lead').
327 %> @param 'UseHamiltonian' Set true if the interface region is matched to the whole Hamiltonian of the scattering center, or false (default) if the surface Green operator of the scattering center is used in the calculations.
328 %> @return [1] The calculated Greens function.
329 %> @return [2] The inverse of the Green operator.
330 %> @return [3] An instance of structure #junction_sites describing the sites in the calculated Green operator.
331  function [Gret, Ginverz, junction_sites] = CustomDysonFunc( obj, varargin ) %NEW output
332 
333  p = inputParser;
334  p.addParameter('gfininv', []);
335  p.addParameter('constant_channels', true);
336  p.addParameter('onlyGinverz', false );
337  p.addParameter('recalculateSurface', [1 2] );
338  p.addParameter('decimate', true );
339  p.addParameter('kulso_szabfokok', []); %The list of sites to be left after the decimation procedure
340  p.addParameter('SelfEnergy', false); %set true to calculate the Dyson equation with the self energy
341  p.addParameter('keep_sites', 'lead'); %Name of sites to be kept (scatter, interface, lead)
342  p.addParameter('UseHamiltonian', false); %true if the interface region is matched to the whole Hamiltonian of the scattering center, false if the surface Green operator of the scattering center is used in the calculations.
343  p.parse(varargin{:});
344  gfininv = p.Results.gfininv;
345  constant_channels = p.Results.constant_channels;
346  onlyGinverz = p.Results.onlyGinverz;
347  recalculateSurface = p.Results.recalculateSurface;
348  decimate = p.Results.decimate;
349  kulso_szabfokok = p.Results.kulso_szabfokok;
350  useSelfEnergy = p.Results.SelfEnergy;
351  keep_sites = p.Results.keep_sites;
352  UseHamiltonian = p.Results.UseHamiltonian;
353 
354 
355  if ~isempty(recalculateSurface)
356 
357  % creating interfaces for the Leads
358  if constant_channels
359  shiftLeads = ones(length(obj.param.Leads),1)*obj.E;
360  else
361  shiftLeads = ones(length(obj.param.Leads),1)*0;
362  end
363 
364  % creating Lead instaces and calculating the retarded surface Green operator/self-energy
365  coordinates_shift = [-2, obj.height+1] + obj.shift;
366  obj.FL_handles.LeadCalc('coordinates_shift', coordinates_shift, 'shiftLeads', shiftLeads, 'transversepotential', obj.transversepotential, ...
367  'SelfEnergy', useSelfEnergy, 'SurfaceGreensFunction', ~useSelfEnergy, 'gauge_field', obj.gauge_field, 'leads', recalculateSurface, 'q', obj.q, ...
368  'leadmodel', obj.leadmodel);
369 
370  for idx = 1:length(recalculateSurface)
371  obj.CreateInterface( recalculateSurface(idx), 'UseHamiltonian', UseHamiltonian );
372  end
373 
374  end
375 
376  [Gret, Ginverz, junction_sites] = CustomDysonFunc@NTerminal(obj, 'gfininv', gfininv, ...
377  'onlyGinverz', onlyGinverz, ...
378  'recalculateSurface', [], ...
379  'decimate', decimate, ...
380  'kulso_szabfokok', kulso_szabfokok, ...
381  'SelfEnergy', useSelfEnergy, ...
382  'keep_sites', keep_sites);
383 
384 
385  end
386 
387 %% CalcFiniteGreensFunction
388 %> @brief Calculates the Green operator of the scattering region by the fast way (see PRB 90, 125428 (2014)).
389 %> @param varargin Cell array of optional parameters identical to #NTerminal.CalcFiniteGreensFunction.
390  function CalcFiniteGreensFunction( obj, varargin )
391 
392  p = inputParser;
393  p.addParameter('gauge_trans', false); % logical: true if want to perform gauge transformation on the Green's function and Hamiltonians
394  p.addParameter('onlyGinv', false);
395  p.parse(varargin{:});
396  gauge_trans = p.Results.gauge_trans;
397  onlyGinv = p.Results.onlyGinv;
398 
399  if ~obj.Scatter_UC.Read('OverlapApplied')
400  obj.Scatter_UC.ApplyOverlapMatrices(obj.E);
401  end
402 
403  % getting the Hamiltonian for the edge slabs
404  K0 = obj.Scatter_UC.Read('K0');
405  K1 = obj.Scatter_UC.Read('K1');
406  K1adj = obj.Scatter_UC.Read('K1adj');
407 
408  q = obj.Scatter_UC.Read('q');
409  if ~isempty( q ) && ~obj.Scatter_UC.Read('HamiltoniansDecimated')
410  K1_transverse = obj.Scatter_UC.Read('K1_transverse');
411  K0 = K0 + K1_transverse*diag(exp(1i*q)) + K1_transverse'*diag(exp(-1i*q));
412  end
413 
414  obj.CreateRibbon();
415 
416  %> transforming the Hamiltonians by SVD if necessary
417  if obj.Scatter_UC.Read('is_SVD_transformed')
418  V = obj.Scatter_UC.Get_V();
419  K0 = V'*K0*V;
420  K1 = V'*K1*V;
421  K1adj = V'*K1adj*V;
422  end
423 
424 
425  %> the first two and last two slabs are added manualy at the end
426  if ~obj.Scatter_UC.Read('is_SVD_transformed')
427  z1 = 1;
428  z2 = obj.height-2;
429  else
430  z1 = 2;
431  z2 = obj.height-3;
432  end
433 
434  obj.G = [];
435  obj.Ginv = [];
436 
437  obj.Scatter_UC.FiniteGreenFunction(z1,z2, 'onlygfininv', true);
438  obj.Ginv = obj.Scatter_UC.Read( 'gfininv' );
439 
440  %> adding to the scattering region the first set of transition layers
441  Neff = obj.Scatter_UC.Get_Neff();
442  non_singular_sites = obj.Scatter_UC.Read( 'kulso_szabfokok' );
443  if isempty( non_singular_sites )
444  non_singular_sites = 1:Neff;
445  end
446 
447  non_singular_sites_edges = [non_singular_sites, size(K0,1)+1:2*size(K0,1)];
448  ginv_edges = -[K0, K1; K1adj, K0];
449  ginv_edges = obj.DecimationFunction( non_singular_sites_edges, ginv_edges );
450 
451  %Ginv = [first_slab, H1, 0;
452  % H1', invG, H1;
453  % 0 , H1', last_slab];
454 
455  obj.Ginv = [ginv_edges(1:Neff,1:Neff), [ginv_edges(1:Neff, Neff+non_singular_sites), zeros(Neff, Neff+size(K0,2))]; ...
456  [ginv_edges(Neff+non_singular_sites, 1:Neff); zeros(Neff, Neff)], obj.Ginv, [zeros(Neff, size(ginv_edges,2)-Neff); ginv_edges(1:Neff, Neff+1:end)]; ...
457  [zeros(size(K0,2),2*Neff), ginv_edges(Neff+1:end, 1:Neff)], ginv_edges(Neff+1:end, Neff+1:end)];
458 
459  [rows, cols] = find( K1 );
460  rows = unique(rows); % cols identical to non_singular_sites
461 
462  %non_singular_sites_Ginv = [1:length(non_singular_sites), size(obj.Ginv,1)-size(K0,1)+1:size(obj.Ginv,1)];
463  non_singular_sites_Ginv = [1:length(non_singular_sites), size(obj.Ginv,1)-size(K0,1)+reshape(rows, 1, length(rows))];
464  obj.Ginv = obj.DecimationFunction( non_singular_sites_Ginv, obj.Ginv );
465 
466 
467  %> terminate the scattering region by the first and last slabs and transform back to the normal space from the SVD representation
468  if obj.Scatter_UC.Read('is_SVD_transformed')
469  %> adding the first and last slab to
470  %Ginv = [first_slab, H1, 0;
471  % H1', invG, H1;
472  % 0 , H1', last_slab];
473 
474  obj.Ginv = [K0, [K1(:, 1:Neff), zeros(size(K0,1), 2*size(K0,2))]; ...
475  [K1adj(1:Neff,:); zeros(size(K0))], obj.Ginv, [zeros(Neff, size(K0,2)); K1]; ...
476  zeros(size(K0)), [zeros(size(K0,1), Neff), K1adj], K0];
477 
478 
479  non_singular_sites_Ginv = [1:size(K0,1), size(obj.Ginv,1)-size(K0,1)+1:size(obj.Ginv,1)];
480  obj.Ginv = obj.DecimationFunction( non_singular_sites_Ginv, obj.Ginv );
481 
482  % transform back to the normal space
483  V_tot = [V, zeros(size(K0)); zeros(size(K0)), V];
484  obj.Ginv = V_tot*obj.Ginv*V_tot';
485  end
486 
487 
488 %disp( obj.Ginv )
489  %> gauge transformation of the vector potential in the effective Hamiltonians
490  if gauge_trans
491  try
492  % gauge transformation on Green's function
493  if ~isempty(obj.Ginv) && ~isempty(obj.PeierlsTransform_Scatter) && isempty(obj.q) && ~isempty(obj.gauge_field)
494  coordinates_scatter = obj.getCoordinates();
495  %> gauge transformation on the inverse Green's function
496  obj.Ginv = obj.PeierlsTransform_Scatter.gaugeTransformation( obj.Ginv, coordinates_scatter, obj.gauge_field );
497  end
498  catch errCause
499  err = MException('EQuUs:Utils:Ribbon:CalcFiniteGreensFunction', 'Unable to perform gauge transformation');
500  err = addCause(err, errCause);
501  save('Error_Ribbon_CalcFiniteGreensFunction.mat')
502  throw(err);
503  end
504  end
505 
506  if ~onlyGinv
507  rcond_Ginv = rcond(obj.Ginv);
508  if isnan(rcond_Ginv ) || abs( rcond_Ginv ) < 1e-15
509  obj.display( 'EQuUs:Utils:Ribbon::CalcFiniteGreensFunction: Regularizing Ginv by SVD', 1);
510  obj.G = obj.inv_SVD( obj.Ginv );
511  else
512  obj.G = inv(obj.Ginv);
513  end
514  obj.Ginv = [];
515  else
516  obj.G = [];
517  end
518 
519  end
520 
521 %% CalcFiniteGreensFunctionFromHamiltonian
522 %> @brief Calculates the Green operator of the scattering region from the whole Hamiltonian.
523 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
524 %> @param 'gauge_trans' Logical value. Set true to perform gauge transformation on the Green operator and on the Hamiltonians.
525 %> @param 'onlyGinv' Logical value. Set true to calculate only the inverse of the surface Greens function #Ginv, or false (default) to calculate #G as well. In the latter case the attribute #Ginv is set to empty at the end.
526 %> @param 'PotInScatter' Obsolete parameter. Use 'scatterPotential' instead.
527 %> @param 'scatterPotential' A function handle pot=f( #coordinates ) or pot=f( #CreateHamiltonians, Energy) for the potential to be applied in the Hamiltonian (used when FiniteGreensFunctionFromHamiltonian=true).
528  function CalcFiniteGreensFunctionFromHamiltonian( obj, varargin )
529 
530  p = inputParser;
531  p.addParameter('gauge_trans', false); % logical: true if want to perform gauge transformation on the Green's function and Hamiltonians
532  p.addParameter('onlyGinv', false);
533  p.addParameter('PotInScatter', []) %OBSOLETE use scatterPotential instead
534  p.addParameter('scatterPotential', []) %NEW overrides optional argument 'PotInScatter', might be a cell array of function handles
535  p.parse(varargin{:});
536  gauge_trans = p.Results.gauge_trans;
537  onlyGinv = p.Results.onlyGinv;
538 
539  scatterPotential = p.Results.PotInScatter;
540  if ~isempty( p.Results.scatterPotential )
541  scatterPotential = p.Results.scatterPotential;
542  end
543 
544  obj.CreateScatter();
545  CreateH = obj.CreateH.CreateClone();
546 
547  Hscatter = CreateH.Read('Hscatter');
548  Hscatter_transverse = CreateH.Read('Hscatter_transverse');
549 
550 
551  %> apply magnetic field in scatter for finite q
552  if ~isempty( obj.PeierlsTransform_Scatter ) && ~isempty(obj.q)
553  obj.display('EQuUs:Utils:Ribbon:CalcFiniteGreensFunctionFromHamiltonian: Applying magnetic field in the scattering region')
554  obj.PeierlsTransform_Scatter.PeierlsTransform( CreateH );
555  end
556 
557  %> apply custom potential in the scattering center
558  if ~isempty(scatterPotential)
559  if iscell( scatterPotential )
560  for idx = 1:length( scatterPotential )
561  obj.ApplyPotentialInScatter( CreateH, scatterPotential{idx} );
562  end
563  else
564  obj.ApplyPotentialInScatter( CreateH, scatterPotential);
565  end
566  Hscatter = CreateH.Read('Hscatter');
567  end
568 
569  non_singular_sites_Ginv = CreateH.Read('kulso_szabfokok');
570 
571 
572  q = CreateH.Read('q');
573  if ~isempty( q ) && ~CreateH.Read('HamiltoniansDecimated')
574  Hscatter = Hscatter + Hscatter_transverse*diag(exp(1i*q)) + Hscatter_transverse'*diag(exp(-1i*q));
575  end
576 
577  obj.display('EQuUs:Utils:Ribbon:CalcFiniteGreensFunctionFromHamiltonian: Calculating the surface Green function of the scattering region.')
578  Hscatter = (sparse(1:size(Hscatter,1), 1:size(Hscatter,1), obj.E, size(Hscatter,1),size(Hscatter,1)) - Hscatter);
579  indexes = false( size(Hscatter,1), 1);
580  indexes(non_singular_sites_Ginv) = true;
581  Hscatter = [ Hscatter( ~indexes, ~indexes) , Hscatter(~indexes, indexes);
582  Hscatter(indexes, ~indexes), Hscatter(indexes, indexes)];
583 
584  obj.G = obj.partialInv( Hscatter, length(non_singular_sites_Ginv) );
585 
586  obj.CreateRibbon('justHamiltonians', true);
587 
588 %disp( inv(obj.G) )
589  %> gauge transformation of the vector potential in the effective Hamiltonians
590  if gauge_trans
591  try
592  % gauge transformation on Green's function
593  if ~isempty(obj.G) && ~isempty(obj.PeierlsTransform_Scatter) && isempty(obj.q) && ~isempty(obj.gauge_field)
594  surface_coordinates = obj.getCoordinates();
595  %> gauge transformation on the inverse Green's function
596  obj.G = obj.PeierlsTransform_Scatter.gaugeTransformation( obj.G, surface_coordinates, obj.gauge_field );
597  end
598  catch errCause
599  err = MException('EQuUs:Utils:Ribbon:CalcFiniteGreensFunctionFromHamiltonian', 'Unable to perform gauge transformation');
600  err = addCause(err, errCause);
601  save('Error_Ribbon_CalcFiniteGreensFunctionFromHamiltonian.mat')
602  throw(err);
603  end
604  end
605 
606  if onlyGinv
607  rcond_G = rcond(obj.G);
608  if isnan(rcond_G ) || abs( rcond_G ) < 1e-15
609  obj.display( 'EQuUs:Utils:Ribbon:CalcFiniteGreensFunctionFromHamiltonian: Regularizing Ginv by SVD',1);
610  obj.Ginv = obj.inv_SVD( obj.G );
611  else
612  obj.Ginv = inv(obj.G);
613  end
614  obj.G = [];
615  else
616  obj.Ginv = [];
617  end
618  end
619 
620 %% CreateRibbon
621 %> @brief Creates the Hamiltonians for the ribbon shaped scattering region.
622 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
623 %> @param 'justHamiltonians' Logical value. Set true to create the Hamiltonian of the unit cell without performing any further calculations. (default value is 'false')
624  function CreateRibbon( obj, varargin )
625 
626  p = inputParser;
627  p.addParameter('justHamiltonians', false);
628  p.parse(varargin{:});
629  justHamiltonians = p.Results.justHamiltonians;
630 
631 
632  if ~obj.Scatter_UC.Read( 'HamiltoniansCreated' )
633  obj.display( 'EQuUs:Utils:Ribbon:CreateRibbon: Creating ribbon Hamiltonian.' )
634  obj.Scatter_UC.CreateHamiltonians( 'toSave', 0);
635  obj.Scatter_UC.ShiftCoordinates( obj.shift );
636 
637  %applying transverse potential
638  obj.ApplyTransversePotential( obj.Scatter_UC )
639  end
640 
641 
642 
643 
644  % apply magnetic field in scatter
645  if ~isempty( obj.PeierlsTransform_Scatter ) && ~obj.Scatter_UC.Read('MagneticFieldApplied') && obj.Opt.magnetic_field_trans_invariant
646  obj.display('EQuUs:Utils:Ribbon:CreateRibbon: Applying magnetic field in ribbon Hamiltonians')
647  obj.PeierlsTransform_Scatter.PeierlsTransformLeads( obj.Scatter_UC );
648  end
649 
650  if justHamiltonians
651  return
652  end
653 
654  if ~obj.Scatter_UC.Read('HamiltoniansDecimated' )
655  obj.display( 'EQuUs:Utils:Ribbon:CreateRibbon: Solving the eigenproblem in the scattering region' )
656 
657  %> Trukkos sajatertek
658  obj.display(['EQuUs:Utils:Ribbon:CreateRibbon: Eigenvalues of the scattering region'])
659  obj.Scatter_UC.TrukkosSajatertekek(obj.E);
660  %> group velocity
661  obj.display(['EQuUs:Utils:Ribbon:CreateRibbon: Group velocity for the scattering region'])
662  obj.Scatter_UC.Group_Velocity();
663  end
664 
665  end
666 
667 
668 %% CreateInterface
669 %> @brief Creates the Hamiltonians for the interface regions between the leads and scattering center.
670 %> @param idx Identification number of the interface region.
671 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
672 %> @param 'UseHamiltonian' Logical value. Set true if the interface region should be created to match to the whole Hamiltonian of the scattering center, false (default) if only the surface Green operator of the scattering center is used in the calculations.
673  function CreateInterface( obj, idx, varargin )
674 
675  p = inputParser;
676  p.addParameter('UseHamiltonian', false); %true if the interface region is matched to the whole Hamiltonian of the scattering center, false if the surface Green operator of the scattering center is used in the calculations.
677  p.parse(varargin{:});
678  UseHamiltonian = p.Results.UseHamiltonian;
679 
680  %> Hamiltoninans of the interface region
681  Interface_Region = obj.Interface_Regions{idx};
682 
683  % The regularization of the interface is performed according to the Leads
684  Leads = obj.FL_handles.Read( 'Leads' );
685  Lead = Leads{idx};
686 
687  Interface_Region.Write( 'coordinates2', obj.getCoordinates() );
688  Interface_Region.Write( 'K0', Lead.Read('K0'));
689  Interface_Region.Write( 'K1', Lead.Read('K1'));
690  Interface_Region.Write( 'K1adj', Lead.Read('K1adj'));
691  Interface_Region.Write( 'K1_transverse', Lead.Read('K1_transverse'));
692  Interface_Region.Write( 'coordinates', Lead.Read('coordinates'));
693  Interface_Region.Write( 'kulso_szabfokok', Lead.Read('kulso_szabfokok'));
694  Interface_Region.Write( 'OverlapApplied', true);
695 
696  coordinates_shift = [1, -1 ]; %relative to the leads
697  Interface_Region.ShiftCoordinates( coordinates_shift(idx) );
698 
699  %> coupling between the interface and the scattering region
700  Surface_sc = obj.createSurface_sc( idx );
701  Surface_sc.ApplyOverlapMatrices(obj.E);
702 
703 
704  Lead_Orientation = Surface_sc.Read('Lead_Orientation');
705  K1 = Surface_sc.Read('K1');
706  K1adj = Surface_sc.Read('K1adj');
707  [rows, cols] = find( K1 );
708  rows = unique(rows);
709  cols = unique(cols); % identical as non_singular_sites
710  if Lead_Orientation == 1
711 
712  if UseHamiltonian
713  Hscatter = obj.CreateH.Read('Hscatter');
714  if isempty( Hscatter)
715  error('EQuUs:utils:Ribbon:CreateInterface: Hamiltonian of the scattering region needs to be constructed first.')
716  end
717  Neff = size(Hscatter,1)-size(K1,2);
718  if ~Lead.Read('is_SVD_transformed')
719  % simple decimation
720  Kcoupling = [K1, sparse([], [], [], size(K1,1), Neff )];
721  Kcouplingadj = [K1adj; sparse([], [], [], Neff, size(K1,1))];
722  else
723  % regularization with SVD
724  Kcoupling = [K1, sparse([], [], [], size(K1,1), Neff )];
725  Kcouplingadj = [K1adj; sparse([], [], [], Neff, size(K1,1))];
726  end
727  else
728  Neff = length(rows); %number of coupling sites at Lead_Oriantation=-1
729  if ~Lead.Read('is_SVD_transformed')
730  % simple decimation
731  Kcoupling = [K1(:,cols), zeros( size(K1,1), Neff )];
732  Kcouplingadj = [K1adj(cols,:); zeros(Neff, size(K1,1))];
733  else
734  % regularization with SVD
735  Kcoupling = [K1(:,cols), zeros(size(K1,1), Neff )];
736  Kcouplingadj = [K1adj(cols,:); zeros(Neff, size(K1,1))];
737  end
738  end
739 
740 
741 
742  elseif Lead_Orientation == -1
743 
744  if UseHamiltonian
745  Hscatter = obj.CreateH.Read('Hscatter');
746  if isempty( Hscatter)
747  error('EQuUs:utils:Ribbon:CreateInterface: Hamiltonian of the scattering region needs to be constructed first.')
748  end
749  Neff = size(Hscatter,1)-size(K1adj,2);
750  if ~Lead.Read('is_SVD_transformed')
751  % simple decimation
752  Kcoupling = [sparse([], [], [], length(cols), Neff ), K1adj(cols,:)];
753  Kcouplingadj = [sparse([], [], [], Neff, length(cols)); K1(:,cols)];
754  else
755  % regularization with SVD
756  Kcoupling = [sparse([], [], [], size(K1adj,1), Neff ), K1adj];
757  Kcouplingadj = [sparse([], [], [], Neff, size(K1,1) ); K1];
758  end
759  else
760  Neff = length(cols); %number of coupling sites at Lead_Oriantation=1
761  if ~Lead.Read('is_SVD_transformed')
762  % simple decimation
763  Kcoupling = [zeros(length(cols), Neff), K1adj(cols,rows)];
764  Kcouplingadj = [zeros(Neff, length(cols)); K1(rows,cols)];
765  else
766  % regularization with SVD
767  Kcoupling = [zeros(size(K1,1), Neff), K1adj(:,rows)];
768  Kcouplingadj = [zeros(Neff, size(K1,2)); K1(rows,:)];
769  end
770  end
771 
772 
773  else
774  error('EQuUs:Utils:Ribbon:CreateInterface', 'Unknown lead orientation');
775  end
776 
777 
778 
779  Interface_Region.Write('Kcoupling', Kcoupling);
780  Interface_Region.Write('Kcouplingadj', Kcouplingadj);
781 
782  % method to adjust the interface region and coupling to the scattering region by an external function.
783  if ~isempty( obj.interfacemodel )
784  obj.interfacemodel( Interface_Region );
785  end
786 
787  Interface_Region.Calc_Effective_Hamiltonians( 0, 'Lead', Lead );
788 
789  end
790 
791 %% CreateClone
792 %> @brief Creates a clone of the present object.
793 %> @return Returns with the cloned object.
794  function ret = CreateClone( obj )
795 
796  ret = Ribbon( 'width', obj.width, ...
797  'height', obj.height, ...
798  'filenameIn', obj.filenameIn, ...
799  'filenameOut', obj.filenameOut, ...
800  'E', obj.E, ...
801  'EF', 0, ...
802  'phi', obj.phi, ...
803  'silent', obj.silent, ...
804  'transversepotential', obj.transversepotential, ...
805  'Opt', obj.Opt, ...
806  'param', obj.param, ...
807  'q', obj.q, ...
808  'leadmodel', obj.leadmodel, ...
809  'interfacemodel', obj.interfacemodel);
810 
811  ret.EF = obj.EF;
812  ret.CreateH = obj.CreateH.CreateClone();
813  ret.FL_handles = obj.FL_handles.CreateClone();
814  ret.Scatter_UC = obj.Scatter_UC.CreateClone();
815  ret.Interface_Regions = cell(size(obj.Interface_Regions));
816  for idx = 1:length(obj.Interface_Regions)
817  ret.Interface_Regions{idx} = obj.Interface_Regions{idx}.CreateClone();
818  end
819  if ~isempty( obj.PeierlsTransform_Scatter )
820  ret.PeierlsTransform_Scatter = obj.PeierlsTransform_Scatter.CreateClone();
821  end
822 
823  if ~isempty( obj.PeierlsTransform_Leads )
824  ret.PeierlsTransform_Leads = obj.PeierlsTransform_Leads.CreateClone();
825  end
826  ret.gauge_field = obj.gauge_field; % function handle for the scalar field to transform the vector potential from Landauy to Landaux
827 
828  end
829 
830 
831 end % methods public
832 
833 methods ( Access = protected )
834 
835 %% setFermiEnergy
836 %> @brief Sets the Fermi energy on the atomic sites for the calculations (use the same units as the elements of the Hamiltonian).
837  function setFermiEnergy( obj )
838  if ~isempty(obj.EF)
839  obj.param.scatter.epsilon = obj.param.scatter.epsilon - obj.EF;
840  for idx = 1:length(obj.param.Leads)
841  obj.param.Leads{idx}.epsilon = obj.param.Leads{idx}.epsilon - obj.EF;
842  end
843  end
844  end
845 
846 %% createSurface_sc
847 %> @brief Creates the copuling Hamiltonians between the scattering and interface region
848 %> @param idx The identification number of the interface region. (Integer value.)
849 %> @return An instance of class #Lead describing the copuling between the scattering and interface region
850  function Surface_sc = createSurface_sc( obj, idx )
851  Surface_sc = obj.Scatter_UC.CreateClone('empty', true);
852 
853  if ~isempty( obj.param.Leads{idx}.vargamma_sc )
854  params = Surface_sc.Read( 'params' );
855  params.vargamma = obj.param.Leads{idx}.vargamma_sc;
856  Surface_sc.Write( 'params', params );
857  end
858 
859  Surface_sc.CreateHamiltonians( 'toSave', 0);
860  if idx == 1
861  coordinates_shift = 0 + obj.shift ;
862  elseif idx == 2
863  coordinates_shift = obj.height-1 + obj.shift;
864  end
865  Surface_sc.ShiftCoordinates( coordinates_shift )
866  Surface_sc.Write('Hanyadik_Lead', idx);
867  Surface_sc.Write('Lead_Orientation', obj.Interface_Regions{idx}.Read('Lead_Orientation'));
868 
869  %> applying transverse potential
870  obj.ApplyTransversePotential( Surface_sc )
871 
872  %> applying magnetic field
873  if ~isempty( obj.PeierlsTransform_Leads )
874  %> In superconducting lead one must not include nonzero magnetic
875  %> field.
876  %> Hamiltonians in transverse computations must remain
877  %> traslational invariant.
878  if ~Surface_sc.isSuperconducting()
879  obj.display('EQuUs:Utils:Ribbon:createSurface_sc: Applying magnetic field in the Hamiltonians')
880  obj.PeierlsTransform_Leads.PeierlsTransformLeads( Surface_sc );
881  else
882  obj.display('EQuUs:Utils:Ribbon:createSurface_sc: Applying gauge transformation in the Hamiltonians')
883  obj.PeierlsTransform_Leads.gaugeTransformationOnLead( Surface_sc, obj.gauge_field );
884  end
885  end
886 
887 
888  end
889 
890 %% ApplyTransversePotential
891 %> @brief Apply the tranvesre potential in the Hamiltonians
892 %> @param Scatter_UC An instance of class #Lead containing the Hamiltonians.
893  function ApplyTransversePotential( obj, Scatter_UC )
894  if ~isempty(obj.transversepotential) && isempty(obj.q) %In transverse computations no transverse potential can be applied
895  coordinates = Scatter_UC.Read('coordinates');
896  if nargin( obj.transversepotential ) == 1
897  potential2apply = obj.transversepotential( coordinates );
898  elseif nargin( obj.transversepotential ) == 2
899  potential2apply = obj.transversepotential( Scatter_UC, obj.E );
900  else
901  error('EQuUs:Utils:Ribbon:ApplyTransversePotential', 'To many input arguments in function handle scatterpotential');
902  end
903 
904  if isprop(coordinates, 'BdG_u')
905  if size( potential2apply, 1) == 1 || size( potential2apply, 2) == 1
906  potential2apply(~coordinates.BdG_u) = -potential2apply(~coordinates.BdG_u);
907  else
908  potential2apply(~coordinates.BdG_u, ~coordinates.BdG_u) = -conj(potential2apply(~coordinates.BdG_u, ~coordinates.BdG_u));
909  end
910  end
911  Scatter_UC.AddPotential( potential2apply );
912  end
913  end
914 
915 
916 
917 
918 %% CreateHandles
919 %> @brief Initializes the attributes of the class.
920  function CreateHandles( obj )
921 
922  CreateHandles@NTerminal( obj )
923 
924  obj.Scatter_UC = obj.FL_handles.SurfaceGreenFunctionCalculator([], 'createCore', 1, 'q', obj.q);
925 
926  end
927 
928 %% calculate_lead_attach_points
929 %> @brief Determines the site indexes at which the leads are connected to the scattering center.
930  function calculate_lead_attach_points( obj )
931  for idx = 1:length(obj.param.Leads)
932  obj.param.Leads{idx}.M = obj.param.scatter.shape.width;
933  end
934  end
935 
936 %% createShape
937 %> @brief Creates the geometry data of the ribbon shaped scattering region.
938  function createShape( obj )
939 
940  if ~isempty( obj.width ) && ~isempty( obj.height )
941  obj.param.scatter.shape.width = obj.width;
942  obj.param.scatter.shape.height = obj.height;
943  end
944 
945  if ~isempty(obj.param.scatter.shape.width)
946  obj.width = obj.param.scatter.shape.width;
947  else
948  err = MException(['EQuUs:utils:', class(obj), ':createShape'], 'Shape is not given correctly, width is missing');
949  throw(err)
950  end
951 
952  if ~isempty(obj.param.scatter.shape.height)
953  obj.height = obj.param.scatter.shape.height;
954  else
955  err = MException(['EQuUs:utils:', class(obj), ':createShape'], 'Shape is not given correctly, height is missing');
956  throw(err)
957  end
958 
959  obj.calculate_lead_attach_points();
960 
961  end
962 
963 end % protected methods
964 
965 
966 methods (Access=protected)
967 
968 
969 %% InputParsing
970 %> @brief Parses the optional parameters for the class constructor.
971 %> @param varargin Cell array of optional parameters (https://www.mathworks.com/help/matlab/ref/varargin.html):
972 %> @param 'width' Integer. The number of the nonsingular atomic sites in the cross section of the ribbon.
973 %> @param 'height' Integer. The height of the ribbon in units of the lattice vector.
974 %> @param 'filenameIn' The input filename containing the computational parameters. (Use parameters 'Op' and 'param' instead)
975 %> @param 'filenameOut' The output filename to export the computational parameters.
976 %> @param 'WorkingDir' The absolute path to the working directoy.
977 %> @param 'E' The energy value used in the calculations (in the same units as the Hamiltonian).
978 %> @param 'EF' The Fermi energy in the same units as the Hamiltonian. Attribute #E is measured from this value. (Use for equilibrium calculations in the zero temperature limit. Overrides the one comming from the external source)
979 %> @param 'silent' Set true to suppress output messages.
980 %> @param 'transversepotential' A function handle pot = f( #coordinates ) or pot=f( #CreateLeadHamiltonians, Energy) of the transverse potential applied in the lead. (Instead of #CreateLeadHamiltonians can be used its derived class)
981 %> @param 'leadmodel' A function handle #Lead=f( idx, E, varargin ) of the alternative lead model with equivalent inputs and return values as #Transport_Interface.SurfaceGreenFunctionCalculator and with E standing for the energy.
982 %> @param 'interfacemodel' A function handle f( #InterfaceRegion ) to manually adjus the interface regions. (Usefull when 'leadmodel' is also given. For example see @InterfaceModel)
983 %> @param 'Opt' An instance of the structure #Opt.
984 %> @param 'param' An instance of the structure #param.
985 %> @param 'q' The transverse momentum quantum number.
986  function InputParsing( obj, varargin )
987 
988  p = inputParser;
989  p.addParameter('width', obj.width);
990  p.addParameter('height', obj.height);
991  p.addParameter('filenameIn', obj.filenameIn, @ischar);
992  p.addParameter('filenameOut', obj.filenameOut, @ischar);
993  p.addParameter('WorkingDir', obj.WorkingDir, @ischar);
994  p.addParameter('E', obj.E, @isscalar);
995  p.addParameter('EF', obj.EF);
996  p.addParameter('silent', obj.silent);
997  p.addParameter('transversepotential', obj.transversepotential);
998  p.addParameter('leadmodel', obj.leadmodel); %individual physical model for the contacts
999  p.addParameter('interfacemodel', obj.interfacemodel); %individual physical model for the interface regions
1000  p.addParameter('Opt', obj.Opt);
1001  p.addParameter('param', obj.param);
1002  p.addParameter('q', obj.q);
1003 
1004  p.parse(varargin{:});
1005 
1006  InputParsing@NTerminal( obj, 'filenameIn', p.Results.filenameIn, ...
1007  'filenameOut', p.Results.filenameOut, ...
1008  'WorkingDir', p.Results.WorkingDir, ...
1009  'E', p.Results.E, ...
1010  'EF', p.Results.EF, ...
1011  'silent', p.Results.silent, ...
1012  'leadmodel', p.Results.leadmodel, ...
1013  'interfacemodel', p.Results.interfacemodel, ...
1014  'Opt', p.Results.Opt, ...
1015  'param', p.Results.param, ...
1016  'q', p.Results.q);
1017 
1018 
1019  obj.width = p.Results.width;
1020  obj.height = p.Results.height;
1021  obj.transversepotential = p.Results.transversepotential;
1022 
1023  end
1024 
1025 end % methdos private
1026 
1027 
1028 
1029 end
A class describing an N-terminal geometry for equilibrium calculations mostly in the zero temperature...
Definition: NTerminal.m:38
function isSuperconducting()
Test, whether the lead is in the superconducting phase or not.
function InputParsing(varargin)
Parses the optional parameters for the class constructor.
lead_param Leads
A list of structures lead_param containing the physical parameters for the scattering region...
Definition: structures.m:49
function Landaux(x, y, eta_B, Aconst, height, lattice_constant)
Vector potential in the Landau gauge parallel to the x direction.
function CreateClone(varargin)
Creates a clone of the present class.
Property params
An instance of the structure lead_param.
Structure Opt contains the basic computational parameters used in EQuUs.
Definition: structures.m:120
Structure open_channels describes the open channel in the individual leads.
Definition: structures.m:219
function Conduktance()
Calculates the conductance matrix from the scattering matrix.
Structure shape contains data about the geometry of the scattering region.
Definition: structures.m:166
function Write(MemberName, input)
Sets the value of an attribute in the class.
Class to create and store Hamiltonian of the translational invariant leads.
A class for calculations on a ribbon of finite width for equilibrium calculations mostly in the zero ...
Definition: Ribbon.m:34
Property H1
The coupling Hamiltonian between the unit cells.
function CalcFiniteGreensFunction(varargin)
Calculates the Green operator of the scattering region.
function Transport(Energy, B)
creating the Ribbon class representing the twoterminal setup
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.
Property varargin
list of optional parameters (see http://www.mathworks.com/help/matlab/ref/varargin.html for details)
function createOutput(filename, Opt, param)
This function creates an output file containing the running parameters.
function Landauy(x, y, eta_B)
Vector potential in the Landau gauge parallel to the y direction.
function CreateHamiltonians(Opt, param, varargin)
Constructor of the class.
function CreateHamiltonians(varargin)
Creates the Hamiltonians H_0 and H_1 of the lead.
Property q
A transverse momentum.
function Reset()
Resets all elements in the class.
Property Hcoupling
Coupling Hamiltonian from the interface region to the scattering region.
A class describing the interface region between the scattering region and a lead. ...
A class to calculate the Green functions and self energies of a translational invariant lead The nota...
Definition: Lead.m:29
Property varargin
list of optional parameters (see http://www.mathworks.com/help/matlab/ref/varargin.html for details)
Structure param contains data structures describing the physical parameters of the scattering center ...
Definition: structures.m:45
Property Lead_Orientation
The orientation of the lead. Set +1 is the "incoming" direction of the propagating states is defined ...
scatter_param scatter
An instance of the structure scatter_param containing the physical parameters for the scattering regi...
Definition: structures.m:47
Structure sites contains data to identify the individual sites in a matrix.
Definition: structures.m:247
function ShiftCoordinates(shift)
Shifts the coordinates of the sites by an integer multiple of the lattice vector Coordinates.a.
function Read(MemberName)
Query for the value of an attribute in the class.
epsilon
A physical parameter, see the individual lattice documentations for details.
Definition: structures.m:87
Property Opt
An instance of structure Opt.
Definition: Messages.m:31
shape shape
An instance of structure shape describing the geometry of the scattering region.
Definition: structures.m:74
function CreateClone(varargin)
Creates a clone of the present Lead object.
function SurfaceGreenFunctionCalculator(idx, varargin)
Calculates the surface Green&#39;s function or the self energy of a Lead.
function gauge_field(x, y, eta_B, Aconst, height, lattice_constant)
Scalar gauge field connecting Landaux and Landauy gauges.
function Conductance2()
Conductance calculated by Eq (19) in PRB 73 085414.
function createShape(theNode)
Property E
The energy value for which the TrukkosSajatertekek eigenvalue problem was solved. ...
width
The number of sites in the cross section.
Definition: structures.m:168
Structure containing the coordinates and other quantum number identifiers of the sites in the Hamilto...
Definition: Coordinates.m:24
function AddPotential(V)
Adds on-site potential to the Hamiltonian H0.
function Read(MemberName)
Query for the value of an attribute in the class.
Structure junction_sites contains data to identify the individual sites of the leads, the interface regions and the scattering region in a matrix.
Definition: structures.m:236
Property coordinates
An instance of the structure coordinates.
A class to create and store Hamiltonian of the scattering region.