2 % Copyright (C) 2009-2015 Peter Rakyta, Ph.D.
4 % This program is free software: you can redistribute it and/or modify
5 % it under the terms of the GNU General Public License as published by
6 % the Free Software Foundation, either version 3 of the License, or
7 % (at your option) any later version.
9 % This program is distributed in the hope that it will be useful,
10 % but WITHOUT ANY WARRANTY; without even the implied warranty of
11 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 % GNU General Public License
for more details.
14 % You should have received a copy of the GNU General Public License
15 % along with
this program. If not, see http:
17 %> @addtogroup basic Basic Functionalities
20 %> @brief Class to create the Hamiltonian of one unit cell in a translational invariant lead made of square lattice structure, including the SSH model.
22 %> @brief Class to create the Hamiltonian of one unit cell in a translational invariant lead made of square lattice structure, including the SSH model.
27 methods (Static =
true)
28 %% SquareLattice_Lead_Hamiltonians
29 %> @brief Creates
Hamiltonians H_0 and H_1
for square lattice structure
31 %> @
param M Number of
sites in the cross section of the lead.
32 %> @
param varargin Optional parameters (https:
33 %> @
param 'q' The transverse momentum. Set to empty (
default)
for computations without transverse momentums.
34 %> @
return [1] The Hamiltonian of one slab in the ribbon.
35 %> @
return [2] The coupling between the slabs.
36 %> @
return [3] The transverse coupling between the slabs
for transverse calculations.
37 %> @
return [4] A structure #coordinates containing the coordinates of the
sites.
38 function [H0, H1, H1_transverse, coordinates] = SquareLattice_Lead_Hamiltonians(
lead_param, M, varargin)
41 error([
'EQuUs:',
class(obj),
':SquareLattice_Lead_Hamiltonians'],
'The input parameter M is empty')
45 p.addParameter(
'q', []);
53 H0 = sparse(1:M,1:M,epsilon,M,M) - sparse(1:M-1,2:M,vargamma,M,M) - sparse(2:M,1:M-1,vargamma,M,M);
54 H1 = sparse(1:M,1:M,-vargamma,M,M);
57 coordinates.x = transpose(1:M);
58 coordinates.y = zeros(M,1);
59 coordinates.a = [0; 1];
71 coordinates.b = [M,0];
72 H1_transverse = sparse(M,1,-vargamma,M,M);
77 %% SSH_Lead_Hamiltonians
78 %> @brief Creates
Hamiltonians H_0 and H_1 of the SSH model.
80 %> @
param varargin Optional parameters (https:
81 %> @
param 'q' The transverse momentum. Set to empty (
default)
for computations without transverse momentums.
82 %> @
return [1] The Hamiltonian of one slab in the ribbon.
83 %> @
return [2] The coupling between the slabs.
84 %> @
return [3] The transverse coupling between the slabs
for transverse calculations.
85 %> @
return [4] A structure #coordinates containing the coordinates of the
sites.
86 function [H0, H1, H1_transverse, coordinates] = SSH_Lead_Hamiltonians(
lead_param, varargin)
89 p.addParameter(
'q', []);
105 H0(1,2) = -vargamma1;
106 H0(2,1) = -vargamma1;
107 H0(1,1) = deltaAB/2 + epsilon;
108 H0(2,2) = -deltaAB/2 + epsilon;
110 H1(2,1) = -vargamma3;
115 coordinates.x = transpose(1:2);
116 coordinates.y = zeros(2,1);
117 coordinates.a = [0; 1];
120 error(
'Transverse computations are not allowed in the SSH model')
126 %% Lieb_Lead_Hamiltonians
127 %> @brief Creates
Hamiltonians H_0 and H_1
for Lieb lattice structure
129 %> @
param M Number of
sites in the cross section of the lead.
130 %> @
param varargin Optional parameters (https:
131 %> @
param 'q' The transverse momentum. Set to empty (
default)
for computations without transverse momentums.
132 %> @
return [1] The Hamiltonian of one slab in the ribbon.
133 %> @
return [2] The coupling between the slabs.
134 %> @
return [3] The transverse coupling between the slabs
for transverse calculations.
135 %> @
return [4] A structure #coordinates containing the coordinates of the
sites.
136 function [H0, H1, H1_transverse, coordinates] = Lieb_Lead_Hamiltonians(
lead_param, M, varargin)
139 error([
'EQuUs:',
class(obj),
':Lieb_Lead_Hamiltonians'],
'The input parameter M is empty')
143 p.addParameter(
'q', []);
145 p.parse(varargin{:});
153 h00=toeplitz([epsilon,-vargamma,0]);
154 h01=-vargamma*[0 0 0; 0 0 0; 0 1 0];
155 h1=-vargamma*[0 1 0; 0 0 0; 0 0 0];
159 H0=sparse([],[],[], N+2,N+2);
160 H1=sparse([],[],[], N+2,N+2);
164 H0 = H0 + sparse(l:3:N,m:3:N,h00(l,m),N+2,N+2) + sparse(l:3:N-3,m+3:3:N,h01(l,m),N+2,N+2) + sparse(l+3:3:N,m:3:N-3,h01(m,l),N+2,N+2);
165 H1 = H1 + sparse(l:3:N,m:3:N,h1(l,m),N+2,N+2);
168 H0(N+1:N+2,N+1:N+2)=[epsilon, -vargamma; -vargamma, epsilon];
171 H1(N+1,N+2)=-vargamma;
174 coordinates.x = zeros(N+2,1);
175 coordinates.x(1:3:N+2) = transpose(1:M+1);
176 coordinates.x(2:3:N+2) = transpose(1:M+1);
177 coordinates.x(3:3:N) = transpose(1:M)+0.5;
178 coordinates.y = zeros(N+2,1);
179 coordinates.y(1:3:N+2) = 0.5;
180 coordinates.a = [0; 1];
185 H0=sparse([],[],[], N,N);
186 H1=sparse([],[],[], N,N);
190 H0 = H0 + sparse(l:3:N,m:3:N,h00(l,m),N,N) + sparse(l:3:N-3,m+3:3:N,h01(l,m),N,N) + sparse(l+3:3:N,m:3:N-3,h01(m,l),N,N);
191 H1 = H1 + sparse(l:3:N,m:3:N,h1(l,m),N,N);
196 coordinates.x = zeros(N,1);
197 coordinates.x(1:3:N) = transpose(1:M);
198 coordinates.x(2:3:N) = transpose(1:M);
199 coordinates.x(3:3:N) = transpose(1:M)+0.5;
200 coordinates.y = zeros(N,1);
201 coordinates.y(1:3:N) = 0.5;
202 coordinates.a = [0; 1];
209 H1_transverse = sparse([],[],[], N,N);
212 H1_transverse = H1_transverse + sparse(3*M-3+l,m,h01(l,m),3*M,3*M);
215 coordinates.b = [M,0];
vargamma3
A physical parameter, see the individual lattice documentations for details.
Class to create the Hamiltonian of one unit cell in a translational invariant lead made of square lat...
function Transport(Energy, B)
creating the Ribbon class representing the twoterminal setup
function Hamiltonians(varargin)
Function to create the custom Hamiltonians for the 1D chain.
Structure containing the physical parameters describing a given lead.
vargamma
A physical parameter, see the individual lattice documentations for details.
Structure param contains data structures describing the physical parameters of the scattering center ...
vargamma1
A physical parameter, see the individual lattice documentations for details.
Structure sites contains data to identify the individual sites in a matrix.
epsilon
A physical parameter, see the individual lattice documentations for details.
deltaAB
A physical parameter, see the individual lattice documentations for details.
function structures(name)