Coverage for kwant/solvers/common.py : 89%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
# Copyright 2011-2014 Kwant authors. # # This file is part of Kwant. It is subject to the license terms in the file # LICENSE.rst found in the top-level directory of this distribution and at # http://kwant-project.org/license. A list of Kwant authors can be found in # the file AUTHORS.rst at the top-level directory of this distribution and at # http://kwant-project.org/authors.
# Until v0.13.0, scipy.sparse did not support making block matrices out of # matrices with one dimension being zero: # https://github.com/scipy/scipy/issues/2127 Additionally, 'scipy.sparse.bmat' # didn't support matrices with zero size until v0.18: # https://github.com/scipy/scipy/issues/5976. For now we use NumPy dense # matrices as a replacement.
# TODO: Once we depend on scipy >= 0.18, code for the special cases can be # removed from _make_linear_sys, _solve_linear_sys and possibly other places # marked by the line "See comment about zero-shaped sparse matrices at the top # of common.py".
"""Solver class for computing physical quantities based on solving a liner system of equations.
`SparseSolver` is an abstract base class. It cannot be instantiated as it does not specify the actual linear solve step. In order to be directly usable, a derived class must implement the methods
- `_factorize` and `_solve_linear_sys`, that solve a linear system of equations
and the following properties:
- `lhsformat`, `rhsformat`: Sparse matrix format of left and right hand sides of the linear system, respectively. Must be one of {'coo', 'csc', 'csr'}.
- `nrhs`: Number of right hand side vectors that should be solved at once in a call to `solve_linear_sys`, when the full solution is computed (i.e. kept_vars covers all entries in the solution). This should be not too big too avoid excessive memory usage, but for some solvers not too small for performance reasons. """
@abc.abstractmethod def _factorized(self, a): """ Return a preprocessed version of a matrix for the use with `solve_linear_sys`.
Parameters ---------- a : a scipy.sparse.coo_matrix sparse matrix.
Returns ------- factorized_a : object factorized lhs to be used with `_solve_linear_sys`. """ pass
@abc.abstractmethod def _solve_linear_sys(self, factorized_a, b, kept_vars): """ Solve the linar system `a x = b`, returning the part of the result indicated in kept_vars.
Parameters ---------- factorized : object The result of calling `_factorized` for the matrix a. b : sparse matrix The right hand side. Its format must match `rhsformat`. kept_vars : slice object or sequence of integers A sequence of numbers of variables to keep in the solution
Returns ------- output : NumPy matrix Solution to the system of equations. """ pass
check_hermiticity=True, realspace=False, *, params=None): """Make a sparse linear system of equations defining a scattering problem.
Parameters ---------- sys : kwant.system.FiniteSystem Low level system, containing the leads and the Hamiltonian of a scattering region. in_leads : sequence of integers Numbers of leads in which current or wave function is injected. energy : number Excitation energy at which to solve the scattering problem. args : tuple, defaults to empty Positional arguments to pass to the ``hamiltonian`` method. Mutually exclusive with 'params'. check_hermiticity : bool Check if Hamiltonian matrices are in fact Hermitian. Enables deduction of missing transmission coefficients. realspace : bool Calculate Green's function between the outermost lead sites, instead of lead modes. This is almost always more computationally expensive and less stable. params : dict, optional Dictionary of parameter names and their values. Mutually exclusive with 'args'.
Returns ------- (lhs, rhs, indices, num_orb) : LinearSys `lhs` is a scipy.sparse.csc_matrix, containing the left hand side of the system of equations. `rhs` is a list of matrices with the right hand side, with each matrix corresponding to one lead mentioned in `in_leads`. `indices` is a list of arrays of variables in the system of equations corresponding to the the outgoing modes in each lead, or the indices of variables, on which a lead defined via self-energy adds the self-energy. Finally `num_orb` is the total number of degrees of freedom in the scattering region.
lead_info : list of objects Contains one entry for each lead. If `realspace=False`, this is an instance of `~kwant.physics.PropagatingModes` with a corresponding format, otherwise the lead self-energy matrix.
Notes ----- All the leads should implement a method `modes` if `realspace=False` and a method `selfenergy`.
The system of equations that is created will be described in detail elsewhere. """
raise ValueError('System contains no leads.') return_norb=True, params=params)[:2]
raise ValueError('System Hamiltonian is not Hermitian. ' 'Use option `check_hermiticity=False` ' 'if this is intentional.')
# Process the leads, generate the eigenvector matrices and lambda # vectors. Then create blocks of the linear system and add them # step by step.
# Construct a matrix of 1's that translates the # inter-cell hopping to a proper hopping # from the system to the lead. for i in interface)]
else u_out.shape[0]) msg = ('Lead {0} has hopping with dimensions ' 'incompatible with its interface dimension.') raise ValueError(msg.format(leadnum))
shape=(iface_orbs.size, lhs.shape[0]))
sp.csc_matrix(np.dot(svd_v, u_out))) else: v_sp = transf vdaguout_sp = transf.T * sp.csc_matrix(u_out) lead_mat = - ulinv_out
format=self.lhsformat)
-np.dot(svd_v, u_in)) else: vdaguin_sp = transf.T * sp.csc_matrix(-u_in)
# defer formation of the real matrix until the proper # system size is known else: else: for i in interface)]
msg = ('Self-energy dimension for lead {0} does not ' 'match the total number of orbitals of the ' 'sites for which it is defined.') raise ValueError(msg.format(leadnum))
lhs.shape) # defer formation of true rhs until the proper system # size is known
# Resize the right-hand sides to be compatible with the full lhs # self-energy lead shape=(lhs.shape[0], l)) else: # lead with modes mats[1].shape[0])
else:
# A lead with no rhs. else: raise RuntimeError('Unknown right-hand side format')
out_leads=None, in_leads=None, check_hermiticity=True, *, params=None): """ Compute the scattering matrix of a system.
An alias exists for this common name: ``kwant.smatrix``.
Parameters ---------- sys : `kwant.system.FiniteSystem` Low level system, containing the leads and the Hamiltonian of a scattering region. energy : number Excitation energy at which to solve the scattering problem. args : tuple, defaults to empty Positional arguments to pass to the ``hamiltonian`` method. Mutually exclusive with 'params'. out_leads : sequence of integers or ``None`` Numbers of leads where current or wave function is extracted. None is interpreted as all leads. Default is ``None`` and means "all leads". in_leads : sequence of integers or ``None`` Numbers of leads in which current or wave function is injected. None is interpreted as all leads. Default is ``None`` and means "all leads". check_hermiticity : ``bool`` Check if the Hamiltonian matrices are Hermitian. Enables deduction of missing transmission coefficients. params : dict, optional Dictionary of parameter names and their values. Mutually exclusive with 'args'.
Returns ------- output : `~kwant.solvers.common.SMatrix` See the notes below and `~kwant.solvers.common.SMatrix` documentation.
Notes ----- This function can be used to calculate the conductance and other transport properties of a system. See the documentation for its output type, `~kwant.solvers.common.SMatrix`.
The returned object contains the scattering matrix elements from the `in_leads` to the `out_leads` as well as information about the lead modes.
Both `in_leads` and `out_leads` must be sorted and may only contain unique entries. """
else: else: raise ValueError("Lead lists must be sorted and " "with unique entries.")
check_hermiticity, False, params=params)
enumerate(linsys.indices) if i in out_leads])
# Do not perform factorization if no calculation is to be done. out_leads, in_leads, check_hermiticity)
# See comment about zero-shaped sparse matrices at the top of common.py. format=self.rhsformat)
out_leads=None, in_leads=None, check_hermiticity=True, *, params=None): """ Compute the retarded Green's function of the system between its leads.
An alias exists for this common name: ``kwant.greens_function``.
Parameters ---------- sys : `kwant.system.FiniteSystem` Low level system, containing the leads and the Hamiltonian of a scattering region. energy : number Excitation energy at which to solve the scattering problem. args : tuple, defaults to empty Positional arguments to pass to the ``hamiltonian`` method. Mutually exclusive with 'params'. out_leads : sequence of integers or ``None`` Numbers of leads where current or wave function is extracted. None is interpreted as all leads. Default is ``None`` and means "all leads". in_leads : sequence of integers or ``None`` Numbers of leads in which current or wave function is injected. None is interpreted as all leads. Default is ``None`` and means "all leads". check_hermiticity : ``bool`` Check if the Hamiltonian matrices are Hermitian. Enables deduction of missing transmission coefficients. params : dict, optional Dictionary of parameter names and their values. Mutually exclusive with 'args'.
Returns ------- output : `~kwant.solvers.common.GreensFunction` See the notes below and `~kwant.solvers.common.GreensFunction` documentation.
Notes ----- This function can be used to calculate the conductance and other transport properties of a system. It is often slower and less stable than the scattering matrix-based calculation executed by `~kwant.solvers.common.smatrix`, and is currently provided mostly for testing purposes and compatibility with RGF code.
It returns an object encapsulating the Green's function elements between the system sites interfacing the leads in `in_leads` and those interfacing the leads in `out_leads`. The returned object also contains a list with self-energies of the leads.
Both `in_leads` and `out_leads` must be sorted and may only contain unique entries. """
else: else: raise ValueError("Lead lists must be sorted and " "with unique entries.") raise ValueError("No output is requested.")
check_hermiticity, True, params=params)
enumerate(linsys.indices) if i in out_leads])
# Do not perform factorization if no calculation is to be done. return GreensFunction(np.zeros((len_kv, len_rhs)), lead_info, out_leads, in_leads, check_hermiticity)
# See comment about zero-shaped sparse matrices at the top of common.py. format=self.rhsformat)
check_hermiticity)
*, params=None): """ Calculate the local density of states of a system at a given energy.
An alias exists for this common name: ``kwant.ldos``.
Parameters ---------- sys : `kwant.system.FiniteSystem` Low level system, containing the leads and the Hamiltonian of the scattering region. energy : number Excitation energy at which to solve the scattering problem. args : tuple of arguments, or empty tuple Positional arguments to pass to the function(s) which evaluate the hamiltonian matrix elements. Mutually exclusive with 'params'. check_hermiticity : ``bool`` Check if the Hamiltonian matrices are Hermitian. params : dict, optional Dictionary of parameter names and their values. Mutually exclusive with 'args'.
Returns ------- ldos : a NumPy array Local density of states at each orbital of the system. """
raise NotImplementedError("ldos for non-Hermitian Hamiltonians " "is not implemented yet.")
# TODO: fix this "self-energy is not implemented yet.")
args, check_hermiticity, params=params)[0]
# Do not perform factorization if no further calculation is needed.
# See comment about zero-shaped sparse matrices at the top of common.py. format=self.rhsformat) slice(linsys.num_orb))
*, params=None): r""" Return a callable object for the computation of the wave function inside the scattering region.
An alias exists for this common name: ``kwant.wave_function``.
Parameters ---------- sys : `kwant.system.FiniteSystem` The low level system for which the wave functions are to be calculated. args : tuple of arguments, or empty tuple Positional arguments to pass to the function(s) which evaluate the hamiltonian matrix elements. Mutually exclusive with 'params'. check_hermiticity : ``bool`` Check if the Hamiltonian matrices are Hermitian. params : dict, optional Dictionary of parameter names and their values. Mutually exclusive with 'args'.
Notes ----- The returned object can be itself called like a function. Given a lead number, it returns a 2d NumPy array that contains the wave function within the scattering region due to each incoming mode of the given lead. Index 0 is the mode number, index 1 is the orbital number.
The modes appear in the same order as the negative velocity modes in `kwant.physics.PropagatingModes`. In Kwant's convention leads are attached so that their translational symmetry points *away* from the scattering region::
left lead SR right lead /---------\ /---\ /---------\ ...-3-2-1-0-X-X-X-0-1-2-3-...
This means that incoming modes (coming from infinity towards the scattering region) have *negative* velocity with respect to the lead's symmetry direction.
Examples -------- >>> wf = kwant.solvers.default.wave_function(some_syst, some_energy) >>> wfs_of_lead_2 = wf(2)
"""
# TODO: figure out what to do with self-energies. ' are not available yet.') args, check_hermiticity, params=params)[0]
slice(self.num_orb))
""" ABC for a linear system solution with variable grouping. """
current_conserving=False):
""" Return slices corresponding to the block from lead_in to lead_out. """
"""Return a slice with the rows in the block corresponding to lead_out. """ self.out_offsets[lead_out + 1])
""" Return a slice with the columns in the block corresponding to lead_in. """ self.in_offsets[lead_in + 1])
"""Return the matrix elements from lead_in to lead_out."""
@abc.abstractmethod def num_propagating(self, lead): """Return the number of propagating modes in the lead.""" pass
@abc.abstractmethod def _transmission(self, lead_out, lead_in): """Return transmission from lead_in to lead_out.
The calculation is expected to be done directly from the corresponding matrix elements. """ pass
"""Return transmission from lead_in to lead_out.
If the option ``current_conserving`` has been enabled for this object, this method will deduce missing transmission values whenever possible.
Current conservation is enabled by default for objects returned by `~kwant.solvers.default.smatrix` and `~kwant.solvers.default.greens_function` whenever the Hamiltonian has been verified to be Hermitian (option ``check_hermiticity``, enabled by default).
""" # OK means: chosen (outgoing, incoming) lead is present.
# Calculate the transmission element by summing over a single # row or column. s_t(*chosen) for chosen[sum_dim] in present[sum_dim])) # Calculate the transmission element at the center of a single # missing row and a single missing column. return (sum(s_t(o, i) - self.num_propagating(o) if o == i else 0 for o, i in product(*present)))
"transmission({0}, {1}).".format(*chosen))
"""Return the conductance matrix.
This is the matrix :math:`C` such that :math:`I = CV` where :math:`I` and :math:`V` are, respectively, the vectors of currents and voltages for each lead.
This matrix is useful for calculating non-local resistances. See Section 2.4 of the book by S. Datta. """ for j in rn] for i in rn]) # Set the diagonal elements such that the sums of columns are zero.
"""A scattering matrix.
Transport properties can be easily accessed using the `~SMatrix.transmission` method (don't be fooled by the name, it can also compute reflection, which is just transmission from one lead back into the same lead.)
`SMatrix` however also allows for a more direct access to the result: The data stored in `SMatrix` is a scattering matrix with respect to lead modes and these modes themselves. The details of this data can be directly accessed through the instance variables `data` and `lead_info`. Subblocks of data corresponding to particular leads are conveniently obtained by `~SMatrix.submatrix`.
`SMatrix` also respects the conservation laws present in the lead, such as spin conservation, if they are declared during system construction. If queried with length-2 sequence the first number is the number of the lead, and the second number is the index of the corresponding conserved quantity. For example ``smatrix.transmission((1, 3), (0, 2))`` is transmission from block 2 of the conserved quantity in lead 0 to the block 3 of the conserved quantity in lead 1.
Attributes ---------- data : NumPy array a matrix containing all the requested matrix elements of the scattering matrix. lead_info : list of data a list containing `kwant.physics.PropagatingModes` for each lead. out_leads, in_leads : sequence of integers indices of the leads where current is extracted (out) or injected (in). Only those are listed for which SMatrix contains the calculated result. """
current_conserving=False): sizes, current_conserving) # in_offsets marks beginnings and ends of blocks in the scattering # matrix corresponding to the in leads. The end of the block # corresponding to lead i is taken as the beginning of the block of # i+1. Same for out leads. For each lead block of the scattering # matrix, we want to pick out the computed blocks of the conservation # law. The offsets of these symmetry blocks are stored in # block_offsets, for all leads. List of lists containing the sizes of # symmetry blocks in all leads. leads_block_sizes[i][j] is the number # of left or right moving modes in symmetry block j of lead # i. len(leads_block_sizes[i]) is the number of blocks in lead i. # If a lead does not have block structure, append None. else: # Symmetry block offsets for all leads - or None if lead does not have # blocks. # Pick out symmetry block offsets for in and out leads np.array(self.block_offsets)[list(self.in_leads)] np.array(self.block_offsets)[list(self.out_leads)] # Block j of in lead i starts at in_block_offsets[i][j]
"""Return a slice with the rows in the block corresponding to lead_out. """ self.out_offsets[lead_out + 1]) else: lead_ind, block_ind = lead_out lead_ind = self.out_leads.index(lead_ind) return slice(self.out_offsets[lead_ind] + self.out_block_offsets[lead_ind][block_ind], self.out_offsets[lead_ind] + self.out_block_offsets[lead_ind][block_ind + 1])
""" Return a slice with the columns in the block corresponding to lead_in. """ self.in_offsets[lead_in + 1]) else: lead_ind, block_ind = lead_in lead_ind = self.in_leads.index(lead_ind) return slice(self.in_offsets[lead_ind] + self.in_block_offsets[lead_ind][block_ind], self.in_offsets[lead_ind] + self.in_block_offsets[lead_ind][block_ind + 1])
"""Return transmission from lead_in to lead_out.
If `lead_in` or `lead_out` are a length-2 sequence, the first number is the number of the lead, and the second number indexes the eigenvalue of the conserved quantity in that lead (e.g. spin) if it was specified. """ # If transmission is between leads and not subblocks, we can use # unitarity (like BlockResult does). else: return self._transmission(lead_out, lead_in)
"""Return the number of propagating modes in the lead.""" return self.sizes[lead]
return ("SMatrix(data=%r, lead_info=%r, out_leads=%r, in_leads=%r)" % (self.data, self.lead_info, self.out_leads, self.in_leads))
""" Retarded Green's function.
Transport properties can be easily accessed using the `~GreensFunction.transmission` method (don't be fooled by the name, it can also compute reflection, which is just transmission from one lead back into the same lead).
`GreensFunction` however also allows for a more direct access to the result: The data stored in `GreensFunction` is the real-space Green's function. The details of this data can be directly accessed through the instance variables `data` and `lead_info`. Subblocks of data corresponding to particular leads are conveniently obtained by `~GreensFunction.submatrix`.
Attributes ---------- data : NumPy array a matrix containing all the requested matrix elements of Green's function. lead_info : list of matrices a list with self-energies of each lead. out_leads, in_leads : sequence of integers indices of the leads where current is extracted (out) or injected (in). Only those are listed for which SMatrix contains the calculated result. """
current_conserving=False): data, lead_info, out_leads, in_leads, sizes, current_conserving)
"""Return the number of propagating modes in the lead.""" self.lead_info[lead].conj().T)
# The number of channels is given by the number of # nonzero eigenvalues of Gamma # rationale behind the threshold from # Golub; van Loan, chapter 5.5.8 eps * np.linalg.norm(gamma, np.inf))
"""Return t * t^dagger in a certain basis."""
# For reflection we have to be more careful self.lead_info[lead_in].conj().T)
# The number of channels is given by the number of # nonzero eigenvalues of Gamma # rationale behind the threshold from # Golub; van Loan, chapter 5.5.8 eps * np.linalg.norm(gamma, np.inf))
return ("GreensFunction(data=%r, lead_info=%r, " "out_leads=%r, in_leads=%r)" % (self.data, self.lead_info, self.out_leads, self.in_leads)) |