Quantum Solvers

The electronic band structure of semiconductor materials is responsible for their light absorption and emission properties as well as for many of their transport properties, ultimately depending on the carriers’ effective masses. These properties are not intrinsic to the material, but depend on external factors, too, most notably the strain and the quantum confinement.

Given the crystalline nature of most semiconductor materials, there will be strain whenever two materials with different crystal lattice constants are grown on top of each other pseudomorphically. Even those with the same lattice constant might be under stress due to other effects such as the presence of impurities or if used at different temperatures having dissimilar thermal expansion coefficients. Quantum confinement, in turn, takes place when the size of the semiconductor material in one or more dimensions is reduced to a few nanometres. In that situation, the energy levels available to the carriers become quantized in the direction of confinement, also changing the density of states. Both conditions take place simultaneously when dealing with strain-balanced quantum wells (QW).

Quantum wells - and more recently quantum wires - have been employed to tune the absorption properties of high efficiency solar cells for the past two decades. The need for appropriate tools to study them in the context of photovoltaics led to the development of the simulation models that were the seed of Solcore. As strained materials with quantum confinement, special care must be taken to obtain a sensible set of parameters for the QW structures, including the band edges with confined energy levels, the effective masses and the absorption coefficient.

Solcore’s approach for evaluating the properties of QWs involves calculating first the effect of strain using a 8-band Pikus-Bir Hamiltonian , treating each material in the structure as bulk, and then using the shifted bands and effective masses to solve the 1D Schödinger equation, after a proper alignment between layers. Finally, the absorption coefficient is calculated based on the 2D density of states, including the effect of excitons.

This section describes the tools used for the calculation of the band structure. The QW absorption profile is discussed in the optics section section.

Bulk 8-band kp solver

solcore.quantum_mechanics.kp_bulk.eight_band_strain_hamiltonian[source]

Hamiltonian to calculate cb, hh, lh, so bands and include strain for the biaxial (along [001]) special case.

See Hamiltonian in ref. but remove row 1 and 6 and column 1 and 6. http://prb.aps.org/pdf/PRB/v73/i12/e125348

solcore.quantum_mechanics.kp_bulk.effective_mass_energy(k, m_eff, E0)[source]
solcore.quantum_mechanics.kp_bulk.kp_bands(material, host_material, kx=0, ky=0, kz=0, return_so=False, graph=False, fit_effective_mass=False, effective_mass_direction='X', additional_traverse=None)[source]
Parameters:
  • material
  • host_material
  • kx
  • ky
  • kz
  • return_so
  • graph
  • fit_effective_mass
  • effective_mass_direction
  • additional_traverse
Returns:

solcore.quantum_mechanics.kp_bulk.KPbands(material, host_material, return_edges_only=False, plot_result=False, t=0, p=1.5707963267948966, fraction=0.2, points=50, vin=None)[source]

New version of the above function that produces either the band edges or the full bands in one specifc direction from Gamma.

  • Calculation of the full Brillouin zone is not allowed. After all, KP is only valid near k=0
  • The calculation of the effective masses is done externally by another function
  • The direction to calculate the bands is provided either with directory angles or a directory vector
  • The number of points of the bands and the fraction up to the border of the Brilluin zone in the given direction are also inputs
  • The default direction is X, as given by the angles theta (t) and phi (p)
solcore.quantum_mechanics.kp_bulk.fit_effective_masses(bands, material, host_material, plot_result=False, dk=0.026)[source]
solcore.quantum_mechanics.kp_bulk.parabolic(k, mr)[source]
solcore.quantum_mechanics.kp_bulk.average_inplane_effective_mass(material, host_material, averaging_points=10, dk=0.026)[source]

Calculates the average in-plane effective mass for the four bands, assuming always the XY plane.

solcore.quantum_mechanics.kp_bulk.kp8x8_bulk(material, host_material, averaging_points=10, dk=0.026)[source]

Combines the above functions to provide a compact result of the band edges and the effective masses.

Quantum well 1D-Schrodinger calculator

solcore.quantum_mechanics.high_level_kp_QW.schrodinger(structure, plot_bands=False, kpoints=40, krange=1000000000.0, num_eigenvalues=10, symmetric=True, quasiconfined=0.0, return_qw_boolean_for_layer=False, Efield=0, blur=False, blurmode='even', mode='kp8x8_bulk', step_size=None, minimum_step_size=0, smallest_feature_steps=20, filter_strength=0, periodic=False, offset=0, graphtype=[], calculate_absorption=False, alpha_params=None, **kwargs)[source]

Solves the Schrodinger equation of a 1 dimensional structure. Depending on the inputs, the method for solving the problem is more or less sophisticated. In all cases, the output includes the band structure and effective masses around k=0 as a function of the position, the energy levels of the QW (electrons and holes) and the wavefunctions, although for the kp4x4 and kp6x6 modes, this is provided as a function of the k value, and therefore there is much more information.

Parameters:
  • structure – The strucutre to solve
  • plot_bands – (False) If the bands should be plotted
  • kpoints
    1. The number of points in the k space
  • krange – (1e-9) The range in the k space
  • num_eigenvalues
    1. Maximum number of eigenvalues to calculate
  • symmetric – (True) If the structure is symmetric, in which case the calculation can be speed up
  • quasiconfined – (0.0 eV) Energy above the band edges that an energy level can have before rejecting it
  • return_qw_boolean_for_layer – (False) Return an boolean array indicating which positions are inside the QW
  • Efield
    1. Electric field.
  • blur – (False) If the potentials and effective masses have to be blurred
  • blurmode – (‘even’) Other values are ‘right’ and ‘left’
  • mode – (‘kp4x4’) The mode of calculating the bands and effective masses. See ‘structure_utilities.structure_to_potentials’
  • step_size – (None) The discretization step of the structure. If none, it is estimated based on the smallest feature.
  • minimum_step_size
    1. The minimum step size.
  • smallest_feature_steps
    1. The number of steps in the smallest feature
  • filter_strength
    1. If > 0, defines the fraction of the wavefunction that has to be inside the QW in order to consider it ‘confined’
  • periodic – (False) If the structure is periodic. Affects the boundary conditions.
  • offset
    1. Energy offset used in the calculation of the energy levels in the case of the ‘bulk’ solvers
  • graphtype – [] If ‘potential’, the band profile and wavefunctions are ploted
Returns:

A dictionary containing the band structure and wavefunctions as a function of the position and k

solcore.quantum_mechanics.structure_utilities.assemble_qw_structure(repeats, well, bulk_l_top, bulk_l_bottom, barrier, well_interlayer=None, structure_label='QW Structure', shift_wells=0)[source]
solcore.quantum_mechanics.structure_utilities.vary_well_widths(structure, fraction=0.5, region_sought='well')[source]

Returns an array of indices which correspond to the region that are quantum wells in the structure.

solcore.quantum_mechanics.structure_utilities.locate_regions(x, structure, region_sought='well')[source]

Returns an array of indices which correspond to the region that are quantum wells in the structure.

solcore.quantum_mechanics.structure_utilities.well_regions(x, structure)[source]
solcore.quantum_mechanics.structure_utilities.text_render(structure, resolution=100)[source]
solcore.quantum_mechanics.structure_utilities.structure_to_potentials(structure, step_size=None, minimum_step_size=0, smallest_feature_steps=20, blur=False, blurmode='even', return_qw_boolean_for_layer=False, Efield=0, mode='kp4x4')[source]

Discretizes the structure as a function of the position, providing the potential for the electrons, HH, LH and SO bands as well as the Luttinger parameters and the effective masses.

The result depend on the chosen mode:
  • kp4x4, calculates the bands and effective masses at k=0 coupling the HH and LH bands
  • kp6x6, calculates the bands and effective masses at k=0 coupling the HH, LH and SO bands
  • kp8x8_bulk, calculates the four bands and fit the effective masses with a parabola around k=0
  • strain, just shifts the bands according to the strain
  • relaxed, do nothing and things are calculated as if we had the bulk, unstrained materials

Notice that although the kp8x8_bulk couples the 4 (8) bands, the effective mass is the result of a fitting around k=0. Therefore, it is not valid to solve 1D problems, but just a bulk-like problem under a parabolic aproximation

Additionlly, the band structure can be blured - to simulate the intermixing of materials - as well as being under an electric field.

solcore.quantum_mechanics.structure_utilities.random() → x in the interval [0, 1).
solcore.quantum_mechanics.potential_utilities.sort_simultaneous(*lists)
solcore.quantum_mechanics.potential_utilities.tridiag_euler(V, z, m, periodic=False, num_eigenvalues=10, quasiconfined=0)[source]

Returns eignvalue and eigenvectors of an arbitrary potential.

A tridiagonal matrix is constructed by writing the variable effective mass Schrodinger equation over a series of mesh points. The eigenvalues of the matrix correspond to the allowed energy levels of the system.

The previous solver, eig, has been replaced by the spare matrix version, eigs, that is faster to compute

solcore.quantum_mechanics.potential_utilities.schroedinger_solve(x, V, m, num_eigenvalues=10, periodic=False, offset=0, electron=True, quasiconfined=0)[source]

Returns normalised wavefuctions from the potential profile.

Arguments:

x – spatial grid V – potential m – effective mass

Keywords:

electron – whether the wavefunctions describe electrons or holes (default: True) num_eigenvalues – Number of eigenvalues to calculate (default = 10) periodic – not to sure what this does (default: False)
solcore.quantum_mechanics.potential_utilities.discard_unconfined_energy(E, psi, V, quasiconfined)[source]
solcore.quantum_mechanics.potential_utilities.discard_unconfined(x, structure, E, psi, threshold=0.8)[source]
solcore.quantum_mechanics.potential_utilities.potentials_to_wavefunctions_energies(x, Ve, me, Vhh, mhh, Vlh, mlh, num_eigenvalues=10, periodic=False, offset=0, filter_strength=0, structure=None, quasiconfined=0, **kwargs)[source]
solcore.quantum_mechanics.potential_utilities.graph(x, Ve, me, Vhh, mhh, Vlh, mlh, **kwargs)[source]
solcore.quantum_mechanics.graphics.L(x, centre, hwhm)[source]
solcore.quantum_mechanics.graphics.structure_graph(structure, **kwargs)[source]
solcore.quantum_mechanics.graphics.normalise_psi(p)[source]
solcore.quantum_mechanics.graphics.prepare_wavefunction_data(x, E, Psi, trim_levels_beyond, linewidth, color, scale, suppress_invert=False, square=False, alpha=1, label='untitled')[source]
solcore.quantum_mechanics.graphics.wavefunctions_graph(x, Ee=None, psi_e=None, Ehh=None, psi_hh=None, Elh=None, psi_lh=None, trim_levels_beyond=0.001, linewidth=1, scale=0.06, colour_by_band=True, **kwargs)[source]
solcore.quantum_mechanics.graphics.potentials_graph(x, Ve=None, Vhh=None, Vlh=None, color='grey', **kwargs)[source]
solcore.quantum_mechanics.graphics.schrodinger_graph(schrodinger_result, **kwargs)[source]
solcore.quantum_mechanics.graphics.split_schrodinger_graph(schrodinger_result, trim_levels_beyond=0.01, linewidth=1, scale=0.03, suppress_invert=False, probability_density=False, wfalpha=0.8, potentialalpha=0.8, **kwargs)[source]
solcore.quantum_mechanics.graphics.split_schrodinger_graph_potentials(schrodinger_result, trim_levels_beyond=0.01, linewidth=1, scale=0.3, suppress_invert=False, probability_density=False, wfalpha=0.8, potentialalpha=0.8, **kwargs)[source]
solcore.quantum_mechanics.graphics.split_schrodinger_graph_LDOS(schrodinger_result, **kwargs)[source]
solcore.quantum_mechanics.graphics.LDOS1D_e(x, E, psi, m, step=0.001, margin=0.02, broad=0.005)[source]
solcore.quantum_mechanics.graphics.LDOS1D_h(x, E, psi, m, step=0.001, margin=0.02, broad=0.005)[source]
solcore.quantum_mechanics.graphics.calculate_in_plane_masses(x, psi, m)[source]

Calculates the in-plane effective mass for each level, considering that the wavefunction leaks into the barriers.

solcore.quantum_mechanics.graphics.prepare_wavefunction_data_only(x, E, Psi, trim_levels_beyond, linewidth, color, scale, suppress_invert=False, square=False, alpha=1, label='untitled')[source]

4-bands kp QW solver (experimental)

solcore.quantum_mechanics.kp_QW.sort_simultaneous(*lists)
solcore.quantum_mechanics.kp_QW.alfa(Eband, m)[source]
solcore.quantum_mechanics.kp_QW.kp4x4(mat, host_material)[source]

Calculates the bandedges and bulk effective masses of a material under biaxial strain at the gamma point. SO band is ignored.

This function provides the position at the zone center of the C, HH and LH bands of a material under biaxial strain. Both, the paralell and transverse effective masses (with respect the Z axis) are provided. SO band is ignored, therefore the masses are independent of strain.

solcore.quantum_mechanics.kp_QW.kp6x6(mat, host_material)[source]

Calculates the bandedges and bulk effective masses of a material under biaxial strain at the gamma point including the SO band.

This function provides the position at the zone center of the C, HH, LH and SO bands of a material under biaxial strain. Both, the paralell and transverse effective masses (with respect the Z axis) are provided.

solcore.quantum_mechanics.kp_QW.fill_hamiltonian_holes_4x4(A1, A2, B1, B2, C, D, delta, block)[source]
solcore.quantum_mechanics.kp_QW.solve_holes_QW_at_kt_4x4(kt, z, fhh, flh, g1, g2, g3, num=(10, 10), quasiconfined=0.0, symmetric=False)[source]
solcore.quantum_mechanics.kp_QW.fill_hamiltonian_holes_6x6(A1, A2, B1, B2, C, D, delta, block)[source]
solcore.quantum_mechanics.kp_QW.solve_holes_QW_at_kt_6x6(kt, z, fhh, flh, g1, g2, g3, num=(10, 10), quasiconfined=0.0, symmetric=False)[source]
solcore.quantum_mechanics.kp_QW.solve_electrons_QW_at_kt_parabolic(kt, z, fe, me, num=(10, 10), quasiconfined=0.0)[source]

Returns eignvalue and eigenvectors of an arbitrary potential.

A tridiagonal matrix is constructed by writing the variable effective mass Schrodinger equation over a series of mesh points. The eigenvalues of the matrix correspond to the allowed energy levels of the system.

The previous solver, eig, has been replaced by the spare matrix version, eigs, that is faster to compute

“Varible effective mass Schordinger equation and tridiagonal solution method.”, “Frensley, W. R. (1991). Numerical evaluation of resonant states. Superlattices and Microstructures, 11(3), 347350. doi:10.1016/0749-6036(92)90396-M”)

solcore.quantum_mechanics.kp_QW.band_sorting(bands, symmetric=False)[source]

Sort the bands in C, HH and LH according to their character at k=0 :param bands: A dictionary containing the bandstructure as calculated by solve_electrons_QW_at_kt_parabolic and

solve_holes_QW_at_kt_4x4
Returns:A dictionary with the same input information but arranged in a different way, including labels
solcore.quantum_mechanics.kp_QW.solve_bandstructure_QW(structure, num=10, kpoints=40, krange=5000000000.0, quasiconfined=0.0, symmetric=False, plot_bands=False)[source]