Eötvös Quantum Utilities  v5.0.144
Providing the Horsepowers in the Quantum Realm
antidot.m
Go to the documentation of this file.
1 %% Eotvos Quantum Transport Utilities - antidot
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 antidot.m
20 %> @brief A class to perform transport calculations on a graphene antidot (i.e., a hollow in a ribbon). Obsolete class, for real calculations a creation of a new class is recommended.
21 %> @}
22 %> @brief A class to perform transport calculations on a graphene antidot (i.e., a hollow in a ribbon). Obsolete class, for real calculations a creation of a new class is recommended.
23 %%
24 classdef antidot < Ribbon
25 
26 properties ( Access = protected )
27  %> A lattice vector in the hexagonal lattice.
28  a1
29  %> A lattice vector in the hexagonal lattice.
30  a2
31  %> Directory name to export the plotted wave functions.
32  waveFncDirnameFull
33  %> The name of the subdirectory to export the plotted wave functions.
34  waveFncSubDirname
35 end
36 
37 properties ( Access = public )
38  %> An instance of structure hole
39  hole
40  %> an instance of structures scatterers
42  %> The strength of the magnetic field in Tesla
43  B
44  %> an instance of structure structures
45  coordinates
46  %> the sites at the edge of the antidot
47  antidot_edge_points
48  %> flux in the hole in units of \phi_0
49  flux_of_hole
50  %> atom-atom distance
51  rCC
52  %> Planck contant
53  h = 6.626e-34;
54  %> The charge of the electron
55  qe = 1.602e-19;
56 end
57 
58 
59 
60 methods ( Access = public )
61 
62 %% antidot
63 %> @brief Constructor of the class.
64 %> @param Opt An instance of the structure #Opt.
65 %> @param interpq An array of transverse momentum points at which fhandle is calculated via interpolation.
66 %> @param varargin Cell array of optional parameters. For details see #InputParsing.
67 %> @return An instance of the class
68  function obj = antidot( varargin )
69  % Racsvektorok
70  obj.a1 = [sqrt(3); -1]*sqrt(3)/2;
71  obj.a2 = [sqrt(3); 1]*sqrt(3)/2;
72 
73  obj.hole = structures('hole');
74  obj.scatterers = structures('scatterers');
75  obj.coordinates = structures('coordinates');
76 
77  obj.B = 10e-10; %in Tesla
78  obj.antidot_edge_points = [];
79  obj.flux_of_hole = [];
80 
81 
82  obj.G = [];
83  obj.Ginv = [];
84 
85  obj.waveFncDirnameFull = [];
86  obj.waveFncSubDirname = [];
87 
88  obj.InputParsing( varargin{:} );
89  obj.generate_geometry();
90 
91  obj.display('Creating antidot object')
92  obj.createShape();
93  obj.setFermiEnergy();
94 
95  createOutput( obj.filenameOut, obj.Opt, obj.param );
96 
97  obj.CreateHandles();
98  % create Hamiltonians and coordinates
99  obj.CreateRibbon();
100  end
101 
102 %% Transport
103 %> @brief Calculates the conductance of an antidot connected to the leads in the "contact/scattering center/contact" arrangement.
104 %> @param Energy The energy value.
105 %> @return C The calculated conductivity
106 %> @return ny Array of the open channel in the leads.
107 %> @return DeltaC Error of the unitarity.
108  function [C,ny,DeltaC] = Transport( obj, Energy )
109  obj.display(' ')
110  obj.display( 'Calculating transport in create_Hole_in_ribbon')
111  obj.setEnergy( Energy );
112  obj.create_scatter_GreensFunction();
113 
114  Dysonfunc = @()obj.CustomDysonFunc('gfininv', obj.Ginv, 'SelfEnergy', false, 'constant_channels', false);
115  try
116  obj.FL_handles.DysonEq( 'CustomDyson', Dysonfunc );
117  catch errCause
118  save('Error_antidot_Transport.mat');
119  for idx = 1:length(errCause.stack)
120  obj.display(errCause.stack(idx), 1)
121  end
122  throw( errCause )
123  end
124  [S,ny] = obj.FL_handles.SmatrixCalc();
125 
126  norma = norm(S*S'-eye(sum(ny)));
127  if norma >= 1e-3
128  obj.display( ['error of the unitarity of S-matrix: ',num2str(norma)] )
129  end
130 
131  if ny(1) ~= ny(2)
132  obj.display( ['openchannels do not match: ',num2str(ny)] )
133  end
134 
135  conductance = obj.FL_handles.Conduktance();
136  conductance = abs([conductance(1,:), conductance(1,:)]);
137  DeltaC = std(conductance);
138 
139  C = mean(conductance);
140 
141  obj.display( ['Energy = ', num2str(Energy), ' conductance = ', num2str(C)])
142 
143 
144  end
145 
146 %% CalcWavefnc
147 %> @brief Calculates the energies and wave functions of the bound states.
148 %> @param Energy The energy value.
149 %> @param closefigure If true, the created figure is closed at the end of the calculations.
150 %> @param varargin Cell array of optional parameters:
151 %> @param 'toPlot' Set true (default) for plotting the wave function, or false otherwise.
152 %> @param 'Einhomegac' Set true if the Energy is given in units of $$\hbar\omega_c$$, or false (default) otherwise.
153 %> @param 'filterHole' Set true (default) for separating antidot bound states from the edge states, or false otherwise.
154 %> @param 'db' The number of the calculated energy eigenvalues. Default is db=30.
155 %> @param 'infinitemass' Set true for using the infinite mass boundary condition, or false (default) otherwise.
156 %> @param 'dotCalc' Set true for a dot, or false (default) otherwise.
157 %> @param 'delta' Every delta-th point of the wave function becomes plotted. Default is delta=4.
158 %> @param 'smoothedge' Set true for using a smooth edge in the calculations, or false (default) otherwise.
159 %> @return ret A structure containing the calculated results.
160 %TODO create structure corresponding to the return value
161  function ret = CalcWavefnc( obj, Energy, closefigure, varargin )
162  obj.display( 'Calculating the wave function' );
163  if ~exist('closefigure', 'var')
164  closefigure = 1;
165  end
166 
167  p = inputParser;
168  p.addParameter('toPlot', 1, @isscalar);
169  p.addParameter('Einhomegac', 0, @isscalar);
170  p.addParameter('filterHole', 1, @isscalar);
171  p.addParameter('db', 30, @isscalar);
172  p.addParameter('infinitemass', 0, @isscalar);
173  p.addParameter('dotCalc', 0, @isscalar);
174  p.addParameter('delta', 4, @isscalar);
175  p.addParameter('smoothedge', 0, @isscalar);
176 
177  p.parse(varargin{:});
178  toPlot = p.Results.toPlot;
179  Einhomegac = p.Results.Einhomegac;
180  filterHole = p.Results.filterHole;
181  db = p.Results.db;
182  infinitemass = p.Results.infinitemass;
183  dotCalc = p.Results.dotCalc;
184  delta = p.Results.delta; %used in wavefunction plotting
185  smoothedge = p.Results.smoothedge;
186 
187  if Einhomegac
188  rCC = obj.rCC*1e-10; %Angstrom
189  phi0 = obj.h/obj.qe; % magnetic flux quantum
190  eta_B = 2*pi/phi0*(rCC)^2*obj.B;
191  homega = 2.97*sqrt(9/2*eta_B);
192  Energy = Energy*homega;
193  end
194 
195 
196 
197  if dotCalc && ~infinitemass
198  filterHole = 0;
199  end
200  [Hscatter, wavefunc_indexes] = obj.create_hole_Hamiltonian( 'infinitemass', infinitemass, 'dotCalc', dotCalc, 'smoothedge', smoothedge );
201  xcoords = obj.coordinates.x';
202  ycoords = obj.coordinates.y';
203 
204  M = obj.width;
205  height = obj.height*2;
206 
207  obj.display(' Calculating eigenvalues and wave functions')
208  [eigvecs,eigvals] = eigs(Hscatter,db,Energy);
209 
210  eigvals = diag( eigvals );
211  [eigvals, IX] = sort(eigvals);
212  eigvecs = eigvecs(:,IX);
213  ret = struct( 'eigenval',[], 'Wavefnc', [], 'xcoords', [], 'ycoords', []);
214 
215  if filterHole
216  filterHoleStates()
217  ret.hole = 1;
218  end
219  obj.display( eigvals)
220  ret.eigenval = eigvals;
221 
222  % reshaping the eigenvectors
223  Wavefnc = cell( size(eigvals ) );
224  for kdx = 1:length(eigvals)
225 
226  eigvec = eigvecs(:,kdx);
227  wavefnc = zeros(size(xcoords));
228 
229  % plotting wavefunction calculated by the faster way
230  sorok = round( (wavefunc_indexes - mod(wavefunc_indexes,M))/M ) + 1;
231  indexes_enrows = logical( mod(wavefunc_indexes,M) == 0 );
232  sorok( indexes_enrows ) = sorok( indexes_enrows ) -1;
233  indexek_tmp = 0;
234  for idx = 1:height
235  indexes_in_row = logical( sorok == idx );
236  indexes_in_row = wavefunc_indexes(indexes_in_row) - M*(idx-1);
237  wavefnc(idx, indexes_in_row ) = eigvec( indexek_tmp + (1:length(indexes_in_row) ) );
238  indexek_tmp = indexek_tmp + length(indexes_in_row);
239  end
240 
241  Wavefnc{kdx} = wavefnc;
242  end
243 
244  ret.eigenval = eigvals;
245  ret.Wavefnc = Wavefnc;
246  ret.xcoords = xcoords;
247  ret.ycoords = ycoords;
248 
249 
250 
251  if ~toPlot
252  return
253  end
254 
255  mkdir( obj.waveFncDirnameFull );
256  mkdir( obj.waveFncDirnameFull, obj.waveFncSubDirname);
257  obj.display( ' Plotting wavefunctions' )
258  for kdx = 1:length(eigvals)
259 
260 
261  figure1 = figure('XVisual',...
262  '0x21 (TrueColor, depth 24, RGB mask 0xff0000 0xff00 0x00ff)',...
263  'Renderer','painters',...
264  'Colormap',[1 1 1;0.948412716388702 0.948412716388702 0.948412716388702;0.896825432777405 0.896825432777405 0.896825432777405;0.845238089561462 0.845238089561462 0.845238089561462;0.793650805950165 0.793650805950165 0.793650805950165;0.742063522338867 0.742063522338867 0.742063522338867;0.690476179122925 0.690476179122925 0.690476179122925;0.638888895511627 0.638888895511627 0.638888895511627;0.58730161190033 0.58730161190033 0.58730161190033;0.581477105617523 0.581477105617523 0.581477105617523;0.575652658939362 0.575652658939362 0.575652658939362;0.569828152656555 0.569828152656555 0.569828152656555;0.564003705978394 0.564003705978394 0.564003705978394;0.558179199695587 0.558179199695587 0.558179199695587;0.552354753017426 0.552354753017426 0.552354753017426;0.546530246734619 0.546530246734619 0.546530246734619;0.540705800056458 0.540705800056458 0.540705800056458;0.534881293773651 0.534881293773651 0.534881293773651;0.52905684709549 0.52905684709549 0.52905684709549;0.523232340812683 0.523232340812683 0.523232340812683;0.517407834529877 0.517407834529877 0.517407834529877;0.511583387851715 0.511583387851715 0.511583387851715;0.505758881568909 0.505758881568909 0.505758881568909;0.499934434890747 0.499934434890747 0.499934434890747;0.494109958410263 0.494109958410263 0.494109958410263;0.488285452127457 0.488285452127457 0.488285452127457;0.482460975646973 0.482460975646973 0.482460975646973;0.476636499166489 0.476636499166489 0.476636499166489;0.470812022686005 0.470812022686005 0.470812022686005;0.464987546205521 0.464987546205521 0.464987546205521;0.459163069725037 0.459163069725037 0.459163069725037;0.445249050855637 0.445249050855637 0.445249050855637;0.431335002183914 0.431335002183914 0.431335002183914;0.417420983314514 0.417420983314514 0.417420983314514;0.403506934642792 0.403506934642792 0.403506934642792;0.389592915773392 0.389592915773392 0.389592915773392;0.375678867101669 0.375678867101669 0.375678867101669;0.361764848232269 0.361764848232269 0.361764848232269;0.347850799560547 0.347850799560547 0.347850799560547;0.333936780691147 0.333936780691147 0.333936780691147;0.320022732019424 0.320022732019424 0.320022732019424;0.306108713150024 0.306108713150024 0.306108713150024;0.292194694280624 0.292194694280624 0.292194694280624;0.278280645608902 0.278280645608902 0.278280645608902;0.264366626739502 0.264366626739502 0.264366626739502;0.25045257806778 0.25045257806778 0.25045257806778;0.236538544297218 0.236538544297218 0.236538544297218;0.222624525427818 0.222624525427818 0.222624525427818;0.208710491657257 0.208710491657257 0.208710491657257;0.194796457886696 0.194796457886696 0.194796457886696;0.180882424116135 0.180882424116135 0.180882424116135;0.166968390345573 0.166968390345573 0.166968390345573;0.153054356575012 0.153054356575012 0.153054356575012;0.139140322804451 0.139140322804451 0.139140322804451;0.12522628903389 0.12522628903389 0.12522628903389;0.111312262713909 0.111312262713909 0.111312262713909;0.0973982289433479 0.0973982289433479 0.0973982289433479;0.0834841951727867 0.0834841951727867 0.0834841951727867;0.0695701614022255 0.0695701614022255 0.0695701614022255;0.0556561313569546 0.0556561313569546 0.0556561313569546;0.0417420975863934 0.0417420975863934 0.0417420975863934;0.0278280656784773 0.0278280656784773 0.0278280656784773;0.0139140328392386 0.0139140328392386 0.0139140328392386;0 0 0]);
265 
266 
267 
268  % Create axes
269  axes1 = axes('Parent',figure1, 'CLim',[0 0.028]);
270  hold on
271 
272  contourf( xcoords(1:delta:end,1:delta:end), ycoords(1:delta:end,1:delta:end), abs(wavefnc(1:delta:end,1:delta:end)), 50, 'LineStyle', 'none','Parent', axes1 );
273  axis equal
274 
275  text(30,30, ['E = ', num2str(eigvals(kdx)*1000),' meV']);
276  set(gcf,'Renderer','painters')
277 
278  filename = ['wavefunc_width',num2str(obj.width),'_height',num2str(obj.height),'_radius',num2str(obj.hole.radius),'_B',num2str(obj.B),'_E',num2str(eigvals(kdx),'%6.5f'),'.eps'];
279  if infinitemass
280  filename = ['infmass_',filename];
281  end
282  if smoothedge
283  filename = ['smoothedge_',filename];
284  end
285  if dotCalc
286  filename = ['dot_', filename];
287  end
288  print('-depsc2', [obj.waveFncDirnameFull, filesep, obj.waveFncSubDirname,'/',filename]);
289 
290 
291  if closefigure
292  close(figure1);
293  end
294 
295 
296 
297  end
298 
299 
300  %----------------------------------------------
301  % nested functions
302 
303  %----------------------------------
304  %this filters out eigenstates related to bound states to the Hole
305  function filterHoleStates()
306 % disp( ' Filtering out hole states' )
307  hole_indexes = true(1,length(eigvals));
308 
309  tol = 1e-2;
310  for iidx = 1:length(hole_indexes)
311  if norm(eigvecs(1:M-20,iidx)) > tol
312  hole_indexes(iidx) = false;
313  end
314  end
315 
316  eigvals = eigvals( hole_indexes );
317  eigvecs = eigvecs(:, hole_indexes );
318  end
319 
320  % end nested functions
321 
322 
323  end
324 
325 %% create_hole_Hamiltonian
326 %> @brief Creates the Hamiltonian for the antidot/dot.
327 %> @param varargin Cell array of optional parameters:
328 %> @param 'infinitemass' Set true for using the infinite mass boundary condition, or false (default) otherwise.
329 %> @param 'dotCalc' Set true for a dot, or false (default) otherwise.
330 %> @param 'smoothedge' Set true for using a smooth edge in the calculations, or false (default) otherwise.
331 %> @return Hscatter The matrix representation of the created Hamiltonian.
332 %> @return wavefunction_indexes The list of the sites included in the antidot/dot system.
333  function [Hscatter, wavefunction_indexes] = create_hole_Hamiltonian( obj, varargin )
334  p = inputParser;
335  p.addParameter('infinitemass', 0, @isscalar);
336  p.addParameter('dotCalc', 0, @isscalar);
337  p.addParameter('smoothedge', 0, @isscalar);
338  p.parse(varargin{:});
339 
340  infinitemass = p.Results.infinitemass;
341  dotCalc = p.Results.dotCalc;
342  smoothedge = p.Results.smoothedge;
343 
345  obj.setEnergy([]);
346 
347  obj.display(' Creating the Hamiltonian of a finite ribbon')
348  obj.CreateScatter( );
349  Hscatter = obj.CreateH.Read( 'Hscatter' );
350 
351  obj.display(' Creating hole in the ribbon ')
352 
353 
354  if infinitemass
355  Hscatter_tmp = Hscatter;
356  V = -2;
357  elseif dotCalc
358  Hscatter_tmp = Hscatter;
359  V = [];
360  else
361  Hscatter_tmp = Hscatter;
362  V = [];
363  end
364 
365 
366  if isempty(obj.antidot_edge_points)
367  obj.getHolePoints();
368  end
369 
370  obj.display( ' cutting out hole points ')
371  M = obj.width;
372  indexes2clear = [];
373  for idx = size(obj.antidot_edge_points,1):-1:1
374  indexek_tmp = (obj.antidot_edge_points(idx)-1)*(2*M) + (obj.antidot_edge_points(idx,2) : obj.antidot_edge_points(idx,3));
375  Hscatter( indexek_tmp,: ) = 0;
376  Hscatter( :, indexek_tmp ) = 0;
377  indexes2clear = [indexes2clear, indexek_tmp];
378  end
379 
380  indexes2clear_tmp = getHoleSitesFromTheOddRows( indexes2clear );
381  indexes2clear = sort( [ indexes2clear, indexes2clear_tmp]);
382 
383 
384  if ~infinitemass
385  indexes_tmp = 1:size(Hscatter_tmp,1);
386  indexes_tmp( indexes2clear) = [];
387 
388  if dotCalc
389  wavefunction_indexes = indexes2clear;
390  Hscatter = Hscatter_tmp;
391  Hscatter( indexes_tmp,: ) = [];
392  Hscatter( :, indexes_tmp ) = [];
393  else
394  if smoothedge
395  wavefunction_indexes = indexes_tmp;
396  Hscatter( indexes2clear,: ) = [];
397  Hscatter( :, indexes2clear ) = [];
398  else
399  wavefunction_indexes = 1:size(Hscatter,1);
400  end
401  end
402  else
403  if dotCalc
404  indexes_tmp = 1:size(Hscatter_tmp,1);
405  indexes_tmp( indexes2clear ) = [];
406  indexes2clear = indexes_tmp;
407  end
408  ApplyInfiniteMass( indexes2clear )
409  wavefunction_indexes = 1:size(Hscatter_tmp,1);
410  end
411 
412  %----------------------------------------------
413  % nested functions
414 
415  %-----------------------------------------
416  function ApplyInfiniteMass( indexek)
417  Hscatter = Hscatter_tmp;
418 
419  sorok = round((indexek - mod(indexek,M))/M)+1;
420  indexes_endrows = logical( mod(indexek,M) == 0 );
421  sorok( indexes_endrows ) = sorok( indexes_endrows ) -1;
422 
423  ps_sorok = logical( mod(sorok,2) == 0 );
424  pt_sorok = ~ps_sorok;
425 
426  ps_sorok = indexek( ps_sorok );
427  pt_sorok = indexek( pt_sorok );
428 
429  idx_A = [ ps_sorok( logical( mod(ps_sorok,2) == 0)), pt_sorok( logical( mod(pt_sorok,2) == 1)) ];
430  idx_B = [ ps_sorok( logical( mod(ps_sorok,2) == 1)), pt_sorok( logical( mod(pt_sorok,2) == 0)) ];
431 
432  Hscatter = Hscatter + sparse( idx_A, idx_A, V, size(Hscatter,1), size(Hscatter,2)) + sparse( idx_B, idx_B, -V, size(Hscatter,1), size(Hscatter,2));
433 
434 
435  end
436 
437 
438  % ----------------------------------------
439  function indexes2clear = getHoleSitesFromTheOddRows( indexek_tmp )
440 
441  indexes2clear = [];
442 
443  for iidx = 1:length( indexek_tmp )
444  related_indexes = find(Hscatter_tmp(indexek_tmp(iidx),:));
445  sorok = round((related_indexes - mod(related_indexes,M))/M)+1;
446  indexes_endrows = logical( mod(related_indexes,M) == 0 );
447  sorok( indexes_endrows ) = sorok( indexes_endrows ) -1;
448 
449  ps_sorok = logical( mod(sorok,2) == 0);
450  related_indexes = related_indexes( ps_sorok );
451 
452  for iiidx = 1:length(related_indexes)
453  idx_tmp = find(Hscatter(related_indexes(iiidx),:));
454  idx_tmp2 = find( idx_tmp ~= related_indexes(iiidx) );
455  if length( idx_tmp2 ) <= 1;
456  if isempty( find(indexes2clear == related_indexes(iiidx),1) )
457  indexes2clear = [indexes2clear, related_indexes(iiidx)];
458  end
459  end
460  end
461 
462  end
463 
464 
465  end
466 
467  % end nested functions
468 
469 
470  end
471 
473 %> @brief Creates function handles of the vector potentials and apply the magnetic filed in the ribbon Hamiltonians.
474  function CreateHandlesForMagneticField( obj )
475  hLandauy = obj.createVectorPotential( obj.B );
476  obj.setHandlesForMagneticField('scatter', hLandauy, 'lead', hLandauy );
477  end
478 
479 end
480 
481 
482 methods ( Access = protected )
483 
485 %> @brief creates handle for vector potential
486  function hLandauy = createVectorPotential( obj, B )
487 
488  rCC = obj.rCC*1e-10; %Angstrom
489  phi0 = obj.h/obj.qe;
490 
491  eta_B = 2*pi/phi0*(rCC)^2*B;
492  Aconst = 0.1; % constant vektor potential in the systems: stabilizes the numerical computation
493 
494  hLandauy = @(x,y)([zeros(size(x,1),1),eta_B*x]);
495 
496 
497 
498 
499  end
500 
501 %% generate_geometry
502 %> @brief creates geometry data of the hole in the ribbon
503  function generate_geometry( obj )
504 
505  obj.display(' ')
506  obj.display(' Generating geometry in antidot ')
507  obj.display( ' ' )
508 
509  % creating geometry with the
510  obj.convert_site_indexes( )
511  end
512 
513 end % methods private
514 
515 methods ( Access = public )
516 
517 %% getHolePoints
518 %> @brief Determines the sites that should be cut off from the ribbon in order to create the hole.
519 %> @returns antidot_edge_points A matrix with culomns z (slab index), zpoints_min, zpoints_max (lower and upper bounds of the sites in the slab z).
520  function antidot_edge_points = getHolePoints( obj )
521 
522  if isempty(obj.Scatter_UC)
523  obj.CreateRibbon('justHamiltonians', true);
524  end
525 
526  M = obj.Scatter_UC.Read( 'M' );
527 
528  deltax_ps = 2; % 2n-1 es 2n kozotti tavolsag rCC egysegekben
529  deltax_pt = 1; % 2n es 2n+1 kozotti tavolsag rCC egysegekben
530  deltax = 3; % 2n-1 es 2n+1 kozotti tavolsag rCC egysegekben
531  deltay = sqrt(3);
532 
533  hole.center = [obj.hole.center(1)*(deltax_ps+deltax_pt)/2; obj.hole.center(2)*deltay];
534 
535  rA2 = obj.a2-obj.a1;
536  hole.radius = obj.hole.radius*(norm(rA2))/2;
537 
538  zmin = round((hole.center(2) - hole.radius)/deltay) - 1;
539  zmax = round((hole.center(2) + hole.radius)/deltay) + 1;
540 
541  antidot_edge_points = zeros( zmax-zmin, 3);
542  hole_edge_idx = 1;
543 
544  for ztmp = zmin:zmax
545  minfound = 0;
546  for idx = 1:M
547  dist = distance_from_center(ztmp, idx );
548  if ( dist<hole.radius && ~minfound )
549  idx_min = idx;
550  minfound = 1;
551  end
552 
553  if ( dist>hole.radius && minfound )
554  idx_max = idx-1;
555  antidot_edge_points(hole_edge_idx,:) = [ztmp, idx_min, idx_max];
556  hole_edge_idx = hole_edge_idx + 1;
557  break;
558  end
559 
560  end
561  end
562 
563  antidot_edge_points(hole_edge_idx:end,:) = [];
564  obj.antidot_edge_points = antidot_edge_points;
565 
566  %----------------------------------------------
567  % nested functions
568 
569  %----------------------------------------------
570  function ret = distance_from_center(z, idx )
571 
572  maradek = mod(idx,2);
573 
574  if maradek==1
575  xcoord = (idx-1)/2*deltax;
576  else
577  xcoord = (idx-2)/2*deltax + deltax_ps;
578  end
579 
580  ycoord = z*deltay;
581  ret = sqrt( (xcoord-hole.center(1))^2 + (ycoord - hole.center(2))^2 );
582 
583  end
584 
585  % end nested functions
586 
587 
588 
589  end
590 
591 %% create_scatter_GreensFunction
592 %> @brief Calculates the surface Greens function of the antidot. Sets the attributes #gfin and #gfininv to the calculated results.
593  function create_scatter_GreensFunction( obj )
594  obj.CreateRibbon();
595  obj.Scatter_UC.ApplyOverlapMatrices(obj.E);
596  K0 = obj.Scatter_UC.Read('K0');
597  K1 = obj.Scatter_UC.Read('K1');
598  K1adj = obj.Scatter_UC.Read('K1adj');
599 
600 
601  % Trukkos sajatertek
602  obj.display(['EQuUs:Utils:', class(obj), ':CalcFiniteGreensFunction: Eigenvalues of the scattering region'])
603  obj.Scatter_UC.TrukkosSajatertekek(obj.E);
604  % group velocity
605  obj.display(['EQuUs:Utils:', class(obj), ':CalcFiniteGreensFunction: Group velocity for the scattering region'])
606  obj.Scatter_UC.Group_Velocity();
607 
608 
609 
610  if isempty(obj.antidot_edge_points)
611  obj.getHolePoints();
612  end
613 
614  %> the first two and last two slabs are added manualy at the end
615  z1 = 1;
616  z2 = obj.height-2;
617 
618 
619  obj.G = GreensFuncAtRibbonEnds();
620  antidot_edge_points = obj.antidot_edge_points;
621  scatterer_points = getScattererPoints();
622 
623  if ~isempty( scatterer_points )
624  ghole = GreensFuncAtHoleEdge();
625  [ghopp1, ghopp2] = GreensFuncAtCouplings();
626  obj.G = [obj.G,ghopp1; ghopp2,ghole];
627  end
628 
629  if isnan(rcond(obj.G) ) || abs( rcond(obj.G) ) < 1e-15
630  obj.Ginv= obj.inv_SVD( obj.G );
631  else
632  obj.Ginv= inv(obj.G);%gfin\eye(size(gfin));%inv(gfin);
633  end
634  obj.G = [];
635 
636  if ~obj.scatterers.aremissing
637  %applying scattering potentials
638  for ldx = 1:length(obj.scatterers.siteindexes)
639  obj.Ginv( obj.scatterers.siteindexes(ldx)+4*M, obj.scatterers.siteindexes(ldx)+4*M) = ...
640  obj.Ginv( obj.scatterers.siteindexes(ldx)+4*M, obj.scatterers.siteindexes(ldx)+4*M) - obj.scatterers.potentials(ldx);
641  end
642  obj.Ginv= obj.Ginv( [1+M:3*M,obj.scatterers.siteindexes+4*M], [1+M:3*M,obj.scatterers.siteindexes+4*M] );
643  kulso_szabfokok = 1:2*M;
644  obj.Ginv= obj.DecimationFunction( kulso_szabfokok, gfininv);
645  else
646  obj.Ginv= obj.Ginv( 1+M:3*M, 1+M:3*M );
647  end
648 
649  % adding to the scattering region the first and last slabs
650  Neff = obj.Scatter_UC.Get_Neff();
651  non_singular_sites = obj.Scatter_UC.Read( 'kulso_szabfokok' );
652  if isempty( non_singular_sites )
653  non_singular_sites = 1:Neff;
654  end
655 
656  % getting the Hamiltonian for the edge slabs
657  non_singular_sites_edges = [non_singular_sites, size(K0,1)+1:2*size(K0,1)];
658  ginv_edges = -[K0, K1; K1adj, K0];
659  ginv_edges = obj.DecimationFunction( non_singular_sites_edges, ginv_edges );
660 
661  %Ginv = [first_slab, H1, 0;
662  % H1', invG, H1;
663  % 0 , H1', last_slab];
664 
665  obj.Ginv = [ginv_edges(1:Neff,1:Neff), [ginv_edges(1:Neff, Neff+non_singular_sites), zeros(Neff, Neff+size(K0,2))]; ...
666  [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)]; ...
667  [zeros(size(K0,2),2*Neff), ginv_edges(Neff+1:end, 1:Neff)], ginv_edges(Neff+1:end, Neff+1:end)];
668 
669  [rows, cols] = find( K1 );
670  rows = unique(rows); % cols identical to non_singular_sites
671 
672  %non_singular_sites_Ginv = [1:length(non_singular_sites), size(obj.Ginv,1)-size(K0,1)+1:size(obj.Ginv,1)];
673  non_singular_sites_Ginv = [1:length(non_singular_sites), size(obj.Ginv,1)-size(K0,1)+reshape(rows, 1, length(rows))];
674  obj.Ginv = obj.DecimationFunction( non_singular_sites_Ginv, obj.Ginv );
675 
676  % apply gauge transformation if given
677  coordinates_scatter = obj.getCoordinates();
678  if ~isempty( obj.gauge_field )
679  try
680  % gauge transformation on Green's function
681  if ~isempty(obj.Ginv) && ~isempty(obj.PeierlsTransform_Scatter) && isempty(obj.q)
682  % gauge transformation on the inverse Green's function
683  obj.Ginv= obj.PeierlsTransform_Scatter.gaugeTransformation( obj.Ginv, coordinates_scatter, obj.gauge_field );
684  end
685  catch errCause
686  err = MException('Ribbon:CalcFiniteGreensFunction', 'Unable to perform gauge transformation');
687  err = addCause(err, errCause);
688  save('Error_Ribbon_CalcFiniteGreensFunction.mat')
689  throw(err);
690  end
691  end
692 
693  %----------------------------------------------
694  % nested functions
695 
696 
697  %--------------------------------------------------------------------------
698  % get sites of hole and scatterer that should be obtained in the Green's function
699  function scatterer_points = getScattererPoints()
700 
701  scatterer_points = cell(size(antidot_edge_points,1),3);
702 
703  if ~isempty(antidot_edge_points)
704  point_num_down = antidot_edge_points(1,3)-antidot_edge_points(1,2) + 1;
705  point_num_up = antidot_edge_points(end,3)-antidot_edge_points(end,2) + 1;
706  end
707 
708  for idx = 1:size( antidot_edge_points,1 )
709  ztmp = antidot_edge_points( idx );
710 
711  if idx == 1
712  zpoints = antidot_edge_points(idx,2):antidot_edge_points(idx,3);
713  indexek = 1:point_num_down;
714  elseif idx == size(antidot_edge_points,1)
715  zpoints = antidot_edge_points(idx,2):antidot_edge_points(idx,3);
716  indexek = max(indexek)+1:max(indexek)+point_num_up;
717  else
718  if antidot_edge_points(idx,2) + 6 < antidot_edge_points(idx,3)
719  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,2)+2,antidot_edge_points(idx,2)+3,antidot_edge_points(idx,3)-3,antidot_edge_points(idx,3)-2,antidot_edge_points(idx,3)-1,antidot_edge_points(idx,3)];
720  indexek = max(indexek)+1:max(indexek)+8;
721  elseif antidot_edge_points(idx,2) + 5 < antidot_edge_points(idx,3)
722  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,2)+2,antidot_edge_points(idx,2)+3,antidot_edge_points(idx,3)-2,antidot_edge_points(idx,3)-1,antidot_edge_points(idx,3)];
723  indexek = max(indexek)+1:max(indexek)+7;
724  elseif antidot_edge_points(idx,2) + 4 < antidot_edge_points(idx,3)
725  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,2)+2,antidot_edge_points(idx,3)-2,antidot_edge_points(idx,3)-1,antidot_edge_points(idx,3)];
726  indexek = max(indexek)+1:max(indexek)+6;
727  elseif antidot_edge_points(idx,2) + 3 < antidot_edge_points(idx,3)
728  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,2)+2,antidot_edge_points(idx,3)-1,antidot_edge_points(idx,3)];
729  indexek = max(indexek)+1:max(indexek)+5;
730  elseif antidot_edge_points(idx,2) + 2 < antidot_edge_points(idx,3)
731  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,3)-1,antidot_edge_points(idx,3)];
732  indexek = max(indexek)+1:max(indexek)+4;
733  elseif antidot_edge_points(idx,2) + 1 < antidot_edge_points(idx,3)
734  zpoints = [antidot_edge_points(idx,2),antidot_edge_points(idx,2)+1,antidot_edge_points(idx,3)];
735  indexek = max(indexek)+1:max(indexek)+3;
736  else
737  zpoints = antidot_edge_points(idx,2:3);
738  indexek = max(indexek)+1:max(indexek)+2;
739  end
740 
741  end
742 
743  scatterer_points{idx,1} = antidot_edge_points(idx,1);
744  scatterer_points{idx,2} = zpoints;
745  scatterer_points{idx,3} = ones(size(zpoints)); %1 for hole edge points, 0 otherwise
746 
747  end
748 
749 
750  z = antidot_edge_points(:,1);
751 
752  for idx = 1:length( obj.scatterers.z )
753  index = find( z == round(obj.scatterers.z(idx)/2),1 );
754  if isempty( index )
755  is_hole_edge = 0; %1 for hole edge points, 0 otherwise
756  scatterer_points = [scatterer_points; {round(obj.scatterers.z(idx)/2), obj.scatterers.zpoints(idx), is_hole_edge} ];
757  z = [z; round(obj.scatterers.z(idx)/2)];
758  else
759  scatterer_points{index,2} = [scatterer_points{index,2}, obj.scatterers.zpoints(idx)];
760  is_hole_edge = 0; %1 for hole edge points, 0 otherwise
761  scatterer_points{index,3} = [scatterer_points{index,3}, is_hole_edge];
762  end
763  end
764 
765  getScattererSiteIndexes();
766 
767 
768 
769  %-----------------------------------------------------------
770  function getScattererSiteIndexes()
771  obj.scatterers.siteindexes = zeros( 1,length( obj.scatterers.z ) );
772 
773  index_sum = 0;
774  kkdx = 0;
775  for iidx = 1:size(scatterer_points,1)
776 
777  scatterer_idx = find( scatterer_points{iidx,3}==0 );
778  if ~isempty( scatterer_idx )
779  obj.scatterers.siteindexes(kkdx+1:kkdx+length(scatterer_idx) ) = index_sum + scatterer_idx;
780  kkdx = kkdx + length(scatterer_idx);
781  end
782  index_sum = index_sum + length(scatterer_points{iidx,3});
783  end
784  end
785 
786  end
787 
788 
789 
790  %-------------------------------------------------
791  % calculate the Greens function for the hole edge points
792  function ghole = GreensFuncAtHoleEdge( )
793 
794  point_num = 0;
795  for idx = 1:size(scatterer_points,1)
796  point_num = point_num + length(scatterer_points{idx,2});
797  end
798 
799  ghole = zeros(point_num,point_num);
800  indexek1 = 0;
801  for idx = 1:size( scatterer_points,1 )
802  z1tmp = scatterer_points{idx,1};
803  z1points = scatterer_points{idx,2};
804  indexek1 = max(indexek1)+1:max(indexek1)+length(z1points);
805 
806  indexek2 = 0;
807  for jdx = 1:size( scatterer_points,1 )
808  z2tmp = scatterer_points{jdx,1};
809  z2points = scatterer_points{jdx,2};
810  indexek2 = max(indexek2)+1:max(indexek2)+length(z2points);
811 
812 
813  gtmp = obj.Scatter_UC.InfGreenFunction(z1tmp, z2tmp, 'z1points', z1points, 'z2points', z2points);
814 
815  ghole(indexek1, indexek2) = gtmp;
816  end
817  end
818 
819  ghole(max(indexek1)+1:end,:) = [];
820  ghole(:,max(indexek2)+1:end) = [];
821 
822  end
823 
824 
825  %-------------------------------------------------
826  % calculate the Greens function between the ribbon and points and hole edge
827  function [ghopp1, ghopp2] = GreensFuncAtCouplings()
828 
829 
830  point_num = size( ghole,1 );
831  M = obj.Scatter_UC.Read( 'M' );
832  ghopp1 = zeros(4*M,point_num);
833  ghopp2 = zeros(point_num,4*M);
834 
835  indexek = 0;
836  for idx = 1:size( scatterer_points,1 )
837  ztmp = scatterer_points{idx,1};
838  zpoints = scatterer_points{idx,2};
839  indexek = max(indexek)+1:max(indexek)+length(zpoints);
840 
841  gtmp = obj.Scatter_UC.InfGreenFunction(z1-1, ztmp, 'z2points', zpoints);
842  ghopp1(1:M, indexek) = gtmp;
843 
844  gtmp = obj.Scatter_UC.InfGreenFunction(z1, ztmp, 'z2points', zpoints);
845  ghopp1(1+M:2*M, indexek) = gtmp;
846 
847  gtmp = obj.Scatter_UC.InfGreenFunction(z2, ztmp, 'z2points', zpoints);
848  ghopp1(1+2*M:3*M, indexek) = gtmp;
849 
850  gtmp = obj.Scatter_UC.InfGreenFunction(z2+1, ztmp, 'z2points', zpoints);
851  ghopp1(1+3*M:4*M, indexek) = gtmp;
852 
853 
854 
855 
856  gtmp = obj.Scatter_UC.InfGreenFunction( ztmp, z1-1, 'z1points', zpoints);
857  ghopp2(indexek,1:M) = gtmp;
858 
859  gtmp = obj.Scatter_UC.InfGreenFunction( ztmp, z1, 'z1points', zpoints);
860  ghopp2(indexek,1+M:2*M) = gtmp;
861 
862  gtmp = obj.Scatter_UC.InfGreenFunction( ztmp, z2, 'z1points', zpoints);
863  ghopp2(indexek,1+2*M:3*M) = gtmp;
864 
865  gtmp = obj.Scatter_UC.InfGreenFunction( ztmp, z2+1, 'z1points', zpoints);
866  ghopp2(indexek,1+3*M:4*M) = gtmp;
867 
868  end
869  end
870 
871 
872  %-------------------------------------------------
873  function gfin = GreensFuncAtRibbonEnds()
874 
875  M = obj.Scatter_UC.Read('M');
876  gfin = zeros(4*M,4*M);
877 
878  gtmp = obj.Scatter_UC.InfGreenFunction(z1-1, z1-1);
879  gfin(1:M,1:M) = gtmp;
880  gfin(1+M:2*M, 1+M:2*M) = gtmp;
881  gfin(1+2*M:3*M, 1+2*M:3*M) = gtmp;
882  gfin(1+3*M:4*M, 1+3*M:4*M) = gtmp;
883 
884  gtmp = obj.Scatter_UC.InfGreenFunction(z1-1, z1);
885  gfin(1:M,1+M:2*M) = gtmp;
886  gfin(1+2*M:3*M,1+3*M:4*M) = gtmp;
887 
888  gtmp = obj.Scatter_UC.InfGreenFunction(z1, z1-1);
889  gfin(1+M:2*M,1:M) = gtmp;
890  gfin(1+3*M:4*M, 1+2*M:3*M) = gtmp;
891 
892  gtmp = obj.Scatter_UC.InfGreenFunction(z1-1, z2);
893  gfin(1:M,1+2*M:3*M) = gtmp;
894 
895  gtmp = obj.Scatter_UC.InfGreenFunction(z2, z1-1);
896  gfin(1+2*M:3*M,1:M) = gtmp;
897 
898  gtmp = obj.Scatter_UC.InfGreenFunction(z1-1, z2+1);
899  gfin(1:M,1+3*M:4*M) = gtmp;
900 
901  gtmp = obj.Scatter_UC.InfGreenFunction(z2+1, z1-1);
902  gfin(1+3*M:4*M,1:M) = gtmp;
903 
904  gtmp = obj.Scatter_UC.InfGreenFunction(z1, z2);
905  gfin(1+M:2*M,1+2*M:3*M) = gtmp;
906 
907  gtmp = obj.Scatter_UC.InfGreenFunction(z2, z1);
908  gfin(1+2*M:3*M,1+M:2*M) = gtmp;
909 
910  gtmp = obj.Scatter_UC.InfGreenFunction(z1, z2+1);
911  gfin(1+M:2*M,1+3*M:4*M) = gtmp;
912 
913  gtmp = obj.Scatter_UC.InfGreenFunction(z2+1, z1);
914  gfin(1+3*M:4*M,1+M:2*M) = gtmp;
915 
916  end
917 
918  % end nested functions
919 
920 
921 
922  end
923 
924 
925 end % methods public
926 
927 methods ( Access = private )
928 
929 %% convert_site_indexes
930 %> @brief Converets site indexes into coordinates
931  function convert_site_indexes( obj )
932 
933  xcoord = zeros( obj.width, obj.height*2 );
934  ycoord = zeros( obj.width, obj.height*2 );
935 
936 
937 
938  rA1 = obj.a1+obj.a2;
939  rA2 = obj.a2-obj.a1;
940 
941  radius = obj.hole.radius*(norm(rA2))/2;
942 
943  obj.flux_of_hole = pi*radius^2*(obj.rCC*1e-10)^2*obj.B/(obj.h/obj.qe);
944 
945  rB = [-1;0];
946 
947 
948  for idx = 1:obj.width
949  for jdx = 1:obj.height*2
950 
951  if mod(idx+jdx,2) == 0 % A atom
952  r = double(idx-1)/2*rA1 + double(jdx-1)/2*rA2;
953  else %B atom
954  r = double(idx)/2*rA1 + double(jdx-1)/2*rA2 + rB;
955  end
956 
957  xcoord(idx,jdx) = r(1);
958  ycoord(idx,jdx) = r(2);
959 
960 
961  end
962  end
963 
964  %center = obj.hole.center(1)/2*rA1 + obj.hole.center(2)/2*rA2;
965  %obj.coordinates.outhole = logical( (xcoord - center(1)).^2 + (ycoord - center(2)).^2 >= radius^2 );
966 
967  obj.coordinates.x = xcoord*obj.rCC;
968  obj.coordinates.y = ycoord*obj.rCC;
969 
970 
971 
972 
973  end
974 
975 end % methods private
976 
977 methods ( Access = protected )
978 
979 %% InputParsing
980 %> @brief Parses the optional parameters for the class constructor.
981 %> @param varargin Cell array of optional parameters:
982 %> @param 'width' Integer. The number of the atomic sites in the cross section of the ribbon.
983 %> @param 'height' Integer. The height of the ribbon in units of the lattice vector.
984 %> @param 'filenameIn' Input filename for the xml input structure.
985 %> @param 'filenameOut' Output filename for the xml input structure.
986 %> @param 'WorkingDir' The absolute path to the working directoy.
987 %> @param 'E' The energy value used in the calculations (in the same units as the Hamiltonian).
988 %> @param 'EF' The Fermi energy in the same units as the Hamiltonian. (overrides the one comming from the external source)
989 %> @param 'silent' Set true for suppress the output messages.
990 %> @param 'transversepotential' A function handle pot=f( #coordinates coords) to calculate the transverse potential in the cross section of the ribbon.
991 %> @param 'leadmodel' A function handle #Lead clead=f( idx, E, varargin ) of the alternative lead model with equivalent inputs and return values as class #Lead and with E standing for the energy.
992 %> @param 'interfacemodel' A function handle f( #InterfaceRegion ) to manually adjus the interface regions. (Usefull when 'leadmodel' is also given.)
993 %> @param 'Opt' An instance of the structure #Opt. (Overrides data in the input file if given)
994 %> @param 'param' An instance of the structure #param. (Overrides data in the input file if given)
995 %> @param 'q' The transverse momentum quantum number.
996 %> @param 'B' The strength of the magnetic field in the units of Tesla.
997 %> @param 'hole' An instance of structure #hole.
998 %> @param 'scatterers' An instance of structure #scatterers.
999 %> @param 'radius' The radius of the hole in units of lattice constant.
1000 %> @param 'rCC' The atomic distance.
1001  function InputParsing(obj, varargin)
1002 
1003  p = inputParser;
1004  p.addParameter('width', obj.width);
1005  p.addParameter('height', obj.height);
1006  p.addParameter('filenameIn', [pwd, filesep, 'Basic_Input_zigzag_leads.xml']);
1007  p.addParameter('filenameOut', [pwd, filesep, 'Basic_Input_running_parameters.xml']);
1008  p.addParameter('WorkingDir', pwd);
1009  p.addParameter('E', obj.E);
1010  p.addParameter('EF', obj.EF);
1011  p.addParameter('silent', obj.silent);
1012  p.addParameter('transversepotential', obj.transversepotential);
1013  p.addParameter('leadmodel', obj.leadmodel); %individual physical model for the contacts
1014  p.addParameter('interfacemodel', obj.interfacemodel); %individual physical model for the interface regions
1015  p.addParameter('Opt', []);
1016  p.addParameter('param', []);
1017  p.addParameter('q', []);
1018 
1019  p.addParameter('B', obj.B);
1020  p.addParameter('hole', obj.hole);
1021  p.addParameter('scatterers', obj.scatterers);
1022  p.addParameter('radius', 0);
1023  p.addParameter('rCC', 1.42);
1024 
1025  p.parse(varargin{:});
1026 
1027  InputParsing@Ribbon( obj, 'width', p.Results.width, ...
1028  'height', p.Results.height, ...
1029  'filenameIn', p.Results.filenameIn, ...
1030  'filenameOut', p.Results.filenameOut, ...
1031  'WorkingDir', p.Results.WorkingDir, ...
1032  'E', p.Results.E, ...
1033  'EF', p.Results.EF, ...
1034  'silent', p.Results.silent, ...
1035  'transversepotential', p.Results.transversepotential, ...
1036  'leadmodel', p.Results.leadmodel, ...
1037  'interfacemodel', p.Results.interfacemodel, ...
1038  'Opt', p.Results.Opt, ...
1039  'param', p.Results.param, ...
1040  'q', p.Results.q);
1041 
1042 
1043  obj.B = p.Results.B;
1044  obj.hole = p.Results.hole;
1045  obj.scatterers = p.Results.scatterers;
1046  obj.rCC = p.Results.rCC;
1047 
1048  if p.Results.radius > 0
1049  obj.hole.radius = p.Results.radius;
1050  end
1051 
1052  obj.waveFncDirnameFull = [pwd, filesep, 'wave_function_full'];
1053  obj.waveFncSubDirname = ['width',num2str(obj.width),'_height',num2str(obj.height),'_radius',num2str(obj.hole.radius),'_B',num2str(obj.B)];
1054 
1055 
1056 
1057 
1058  end
1059 
1060 end % methos protected
1061 
1062 
1063 end
1064 
Structure hole contains the data about the antidot used in class antidot.
Definition: structures.m:23
siteindexes
A vector of the site indexes.
Definition: structures.m:123
Structure Opt contains the basic computational parameters used in EQuUs.
Definition: structures.m:60
A class for calculations on a ribbon of finite width for equilibrium calculations mostly in the zero ...
Definition: Ribbon.m:34
function CreateHandlesForMagneticField(B_loc, cRibbon_loc)
Creates and set function handles of the magnetic vector potentials in the Ribbon class.
function Transport(Energy, B)
Calculates the conductance at a given energy value.
function CreateRibbon()
gauge transformation on the inverse Green's function
function Hamiltonians(varargin)
Function to create the custom Hamiltonians for the 1D chain.
A class to perform transport calculations on a graphene antidot (i.e., a hollow in a ribbon)....
Definition: antidot.m:24
function createOutput(filename, Opt, param)
This function creates an output file containing the running parameters.
Structure scatterers contains data on the scattering impurities used in class antidot.
Definition: structures.m:115
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 a1
A lattice vector in the hexagonal lattice.
Definition: antidot.m:31
function CreateHandlesForMagneticField()
Creates function handles of the vector potentials and apply the magnetic filed in the ribbon Hamilton...
Structure param contains data structures describing the physical parameters of the scattering center ...
Definition: structures.m:45
Structure sites contains data to identify the individual sites in a matrix.
Definition: structures.m:187
function gauge_field(x, y, eta_B, Aconst, height, lattice_constant)
Scalar gauge field connecting Landaux and Landauy gauges.
Property E
The energy value for which the TrukkosSajatertekek eigenvalue problem was solved.
function createVectorPotential(B_loc)
Creates the function handle of the magnetic vector potential.
function structures(name)