Core API¶

Core classes for potential energy surfaces, capture coefficient calculations, and quantum mechanics.

Potential¶

The Potential class represents 1D potential energy surfaces with fitting and quantum solving capabilities.

`carriercapture.core.potential.Potential` ¶

Represents a 1D potential energy surface with fitting and quantum solving capabilities.

This class stores potential energy surface data, fits analytical or interpolated functions to the data, and solves the 1D Schrödinger equation to obtain phonon eigenvalues and wavefunctions.

Attributes:

Name	Type	Description
`name`	`str`	Identifier for the potential
`Q_data`	`NDArray[float64] \| None`	Configuration coordinates (sample points), amu^0.5·Å
`E_data`	`NDArray[float64] \| None`	Energy values at sample points, eV
`Q0`	`float`	Minimum configuration coordinate, amu^0.5·Å
`E0`	`float`	Minimum energy, eV
`Q`	`NDArray[float64] \| None`	Interpolated Q grid for evaluation
`E`	`NDArray[float64] \| None`	Energy on Q grid, eV
`fit_type`	`FitType \| None`	Type of fitting function
`fit_func`	`Callable \| None`	Fitted function (Q -> E)
`fit_params`	`Dict[str, Any]`	Hyperparameters for fitting
`eigenvalues`	`NDArray[float64] \| None`	Phonon eigenvalues (eV)
`eigenvectors`	`NDArray[float64] \| None`	Phonon wavefunctions, shape (nev, len(Q))
`nev`	`int`	Number of eigenvalues to compute
`temperature`	`float`	Temperature for thermal weighting (K)
`use_thermal_weights`	`bool`	Enable Boltzmann weighting in fitting

Examples:

Create from file and fit:

>>> pot = Potential.from_file("data.dat", name="DX-center")
>>> pot.fit(fit_type="spline", order=4, smoothness=0.001)
>>> pot.solve(nev=60)
>>> pot.eigenvalues[:5]  # First 5 eigenvalues
array([0.523, 0.569, 0.615, 0.661, 0.707])

Create harmonic potential:

>>> pot = Potential.from_harmonic(hw=0.02, Q0=0.0, E0=0.0)
>>> pot.solve(nev=10)
>>> pot.eigenvalues[:3]  # E_n = ℏω(n + 1/2)
array([0.01, 0.03, 0.05])

Evaluate potential:

>>> pot(5.0)  # Energy at Q=5.0
1.234

`init(name='', Q0=0.0, E0=0.0, nev=30, temperature=293.15, Q_data=None, E_data=None)` ¶

Initialize a Potential.

Parameters:

Name	Type	Description	Default
`name`	`str`	Identifier for the potential	`""`
`Q0`	`float`	Minimum configuration coordinate (amu^0.5·Å)	`0.0`
`E0`	`float`	Minimum energy (eV)	`0.0`
`nev`	`int`	Number of eigenvalues to compute	`30`
`temperature`	`float`	Temperature for thermal weighting (K)	`293.15`
`Q_data`	`NDArray[float64]`	Configuration coordinate data points (amu^0.5·Å)	`None`
`E_data`	`NDArray[float64]`	Potential energy data points (eV)	`None`

`from_file(filepath, name='', resolution=3001, **kwargs)` `classmethod` ¶

Load potential from two-column data file.

File format: Two columns (Q, E) separated by whitespace or comma. Lines starting with '#' are treated as comments.

Parameters:

Name	Type	Description	Default
`filepath`	`str or Path`	Path to data file	required
`name`	`str`	Name for the potential	`""`
`resolution`	`int`	Number of points for interpolation grid	`3001`
`**kwargs`		Additional arguments passed to init	`{}`

Returns:

Name	Type	Description
`pot`	`Potential`	Potential loaded from file

Examples:

>>> pot = Potential.from_file("potential.csv", name="DX-center")
>>> pot.Q_data.shape
(25,)

`from_harmonic(hw, Q0=0.0, E0=0.0, Q_range=(-20.0, 20.0), npoints=5000, **kwargs)` `classmethod` ¶

Create harmonic potential: E = E0 + (1/2) * k * (Q - Q0)².

Where k = (amu/2) * (ℏω/ℏc)²

Parameters:

Name	Type	Description	Default
`hw`	`float`	Phonon energy ℏω (eV)	required
`Q0`	`float`	Equilibrium position (amu^0.5·Å)	`0.0`
`E0`	`float`	Minimum energy (eV)	`0.0`
`Q_range`	`tuple[float, float]`	Range for Q grid (amu^0.5·Å)	`(-20.0, 20.0)`
`npoints`	`int`	Number of grid points	`5000`
`**kwargs`		Additional arguments passed to init	`{}`

Returns:

Name	Type	Description
`pot`	`Potential`	Harmonic potential

Examples:

>>> pot = Potential.from_harmonic(hw=0.02)  # 20 meV phonon
>>> pot.solve(nev=10)
>>> pot.eigenvalues[0]  # Ground state E_0 = ℏω/2
0.01

`fit(fit_type='spline', **fit_kwargs)` ¶

Fit function to data points.

Parameters:

Name	Type	Description	Default
`fit_type`	`FitType`	Fitting method: "spline", "bspline", "harmonic", "polyfunc", "morse", "morse_poly"	`"spline"`
`**fit_kwargs`		Fitting hyperparameters specific to each method: For "spline": - order : int, default=2 Spline degree (1-5) - smoothness : float, default=0.0 Smoothing parameter (0 = interpolating) - weights : array-like, optional Weights for each data point For "harmonic": - hw : float Phonon energy ℏω (eV)	`{}`

Raises:

Type	Description
`ValueError`	If Q_data and E_data are not set
`NotImplementedError`	For fit types not yet implemented

Examples:

>>> pot.fit(fit_type="spline", order=4, smoothness=0.001)
>>> pot.fit(fit_type="harmonic", hw=0.02)

`filter_thermally_accessible(thermal_energy=None)` ¶

Filter data points beyond thermal energy barrier.

Keeps only the "island" of connected data points below the energy threshold that contains the minimum energy point. This is useful for removing high-energy data points that are thermally inaccessible.

Parameters:

Name	Type	Description	Default
`thermal_energy`	`float`	Thermal energy threshold (eV). If None, uses k_B * self.temperature.	`None`

Raises:

Type	Description
`ValueError`	If Q_data or E_data is not set
`ValueError`	If multiple minima are found

Examples:

>>> pot.filter_thermally_accessible(thermal_energy=0.5)  # Keep points within 0.5 eV of minimum

`solve(nev=None, maxiter=None)` ¶

Solve 1D Schrödinger equation for this potential.

Updates self.eigenvalues and self.eigenvectors in place.

Parameters:

Name	Type	Description	Default
`nev`	`int`	Number of eigenvalues to compute. If None, uses self.nev.	`None`
`maxiter`	`int`	Maximum ARPACK iterations. If None, uses default.	`None`

Raises:

Type	Description
`ValueError`	If potential function is not fitted
`RuntimeError`	If ARPACK fails to converge

Examples:

>>> pot.fit(fit_type="spline")
>>> pot.solve(nev=60)
>>> len(pot.eigenvalues)
60

`call(Q)` ¶

Evaluate potential at Q.

Parameters:

Name	Type	Description	Default
`Q`	`float or NDArray[float64]`	Configuration coordinate(s)	required

Returns:

Name	Type	Description
`E`	`float or NDArray[float64]`	Energy at Q

Raises:

Type	Description
`ValueError`	If potential is not fitted

`copy()` ¶

Create a deep copy of this potential.

Returns:

Name	Type	Description
`pot_copy`	`Potential`	Independent copy

`to_dict()` ¶

Serialize potential to dictionary.

Returns:

Name	Type	Description
`data`	`dict`	Dictionary representation

`from_dict(data)` `classmethod` ¶

Deserialize potential from dictionary.

Parameters:

Name	Type	Description	Default
`data`	`dict`	Dictionary representation	required

Returns:

Name	Type	Description
`pot`	`Potential`	Reconstructed potential

ConfigCoordinate¶

The ConfigCoordinate class calculates carrier capture coefficients using configuration coordinate diagrams.

`carriercapture.core.config_coord.ConfigCoordinate` ¶

Configuration coordinate for carrier capture calculations.

Represents the initial and final potential energy surfaces with electron-phonon coupling for capture coefficient calculations.

This class implements the multiphonon carrier capture theory from: Alkauskas et al., Phys. Rev. B 90, 075202 (2014)

Attributes:

Name	Type	Description
`name`	`str`	Identifier
`pot_i`	`Potential`	Initial state potential
`pot_f`	`Potential`	Final state potential
`W`	`float`	Electron-phonon coupling matrix element (eV)
`degeneracy`	`int`	Degeneracy factor
`overlap_matrix`	`NDArray[float64] \| None`	Wavefunction overlaps ⟨ψ_i\|(Q-Q₀)\|ψ_f⟩, shape (nev_i, nev_f)
`delta_matrix`	`NDArray[float64] \| None`	Gaussian energy conservation δ(ΔE), shape (nev_i, nev_f)
`temperature`	`NDArray[float64] \| None`	Temperature grid (K)
`capture_coefficient`	`NDArray[float64] \| None`	Capture coefficient vs temperature (cm³/s)
`partial_capture_coefficient`	`NDArray[float64] \| None`	Detailed capture contributions, shape (nev_i, nev_f, n_temp)

Examples:

Basic usage:

>>> pot_i = Potential.from_harmonic(hw=0.03, Q0=0.0, E0=1.5)
>>> pot_f = Potential.from_harmonic(hw=0.02, Q0=10.0, E0=0.0)
>>> pot_i.solve(nev=180)
>>> pot_f.solve(nev=60)
>>> cc = ConfigCoordinate(pot_i, pot_f, W=0.205, degeneracy=1)
>>> cc.calculate_overlap(Q0=10.0, cutoff=0.25, sigma=0.0075)
>>> cc.calculate_capture_coefficient(volume=1e-21, temperature=np.linspace(100, 500, 50))
>>> cc.capture_coefficient  # cm³/s at each temperature

`init(pot_i, pot_f, name='', W=0.0, degeneracy=1)` ¶

Initialize a ConfigCoordinate.

Parameters:

Name	Type	Description	Default
`pot_i`	`Potential`	Initial state potential	required
`pot_f`	`Potential`	Final state potential	required
`name`	`str`	Identifier for this configuration coordinate	`""`
`W`	`float`	Electron-phonon coupling matrix element (eV)	`0.0`
`degeneracy`	`int`	Degeneracy factor (number of degenerate states)	`1`

`calculate_overlap(Q0, cutoff=0.25, sigma=0.025)` ¶

Calculate phonon wavefunction overlap integrals.

Computes the matrix elements: S_ij = ⟨ψ_i|(Q - Q₀)|ψ_j⟩ = ∫ ψ_i(Q) · (Q - Q₀) · ψ_j(Q) dQ

And the energy-conserving delta function (Gaussian smearing): δ_ij = exp(-ΔE²/(2σ²)) / (σ√(2π))

Only computes overlaps where |E_i - E_j| < cutoff (energy conservation).

Uses vectorization to compute all j overlaps for each i simultaneously, achieving 3-5x speedup over naive nested loops.

Parameters:

Name	Type	Description	Default
`Q0`	`float`	Shift for the coordinate operator (Q - Q₀), typically the equilibrium position used in e-ph coupling calculation (amu^0.5·Å)	required
`cutoff`	`float`	Energy cutoff for filtering overlaps (eV) Overlaps with \|E_i - E_j\| > cutoff are set to zero	`0.25`
`sigma`	`float`	Width of Gaussian delta function (eV) Represents uncertainty in energy conservation	`0.025`

Raises:

Type	Description
`ValueError`	If potentials don't have the same Q grid
`ValueError`	If potentials haven't been solved (no eigenvalues/eigenvectors)

Notes

The overlap integral is computed using trapezoidal rule numerical integration. The energy cutoff significantly reduces computation by skipping pairs with large energy differences (typically ~50% of pairs).

Examples:

>>> cc.calculate_overlap(Q0=10.0, cutoff=0.25, sigma=0.0075)
>>> cc.overlap_matrix.shape
(180, 60)

`calculate_capture_coefficient(volume, temperature)` ¶

Calculate temperature-dependent capture coefficient.

Implements the multiphonon capture rate formula: C(T) = (V·2π/ℏ)·g·W² · Σ_ij occ_i(T) · |S_ij|² · δ_ij

Where: - V: supercell volume - g: degeneracy - W: electron-phonon coupling - occ_i(T): Boltzmann occupation of initial state i - S_ij: overlap integral - δ_ij: energy-conserving delta function

Fully vectorized over temperature axis for efficiency.

Parameters:

Name	Type	Description	Default
`volume`	`float`	Supercell volume where W was calculated (cm³)	required
`temperature`	`NDArray[float64]`	Temperature grid (K), shape (n_temp,)	required

Raises:

Type	Description
`ValueError`	If overlaps haven't been calculated
`ValueError`	If partition function isn't converged (need more eigenvalues)

Notes

The partition function Z(T) = Σ_i exp(-E_i/k_B·T) must be converged, which requires that the occupation of the highest eigenvalue is small: occ(E_max) < 1e-5

If this criterion fails, increase nev for the initial potential.

Examples:

>>> temps = np.linspace(100, 500, 50)
>>> cc.calculate_capture_coefficient(volume=1e-21, temperature=temps)
>>> cc.capture_coefficient  # cm³/s
array([1.23e-10, 2.45e-10, ...])

`to_dict()` ¶

Serialize configuration coordinate to dictionary.

Returns:

Name	Type	Description
`data`	`dict`	Dictionary representation

`from_dict(data)` `classmethod` ¶

Deserialize configuration coordinate from dictionary.

Parameters:

Name	Type	Description	Default
`data`	`dict`	Dictionary representation	required

Returns:

Name	Type	Description
`cc`	`ConfigCoordinate`	Reconstructed configuration coordinate

TransferCoordinate¶

Experimental

The TransferCoordinate class implements Marcus theory for charge transfer. This is an experimental feature (Phase 3) and may have incomplete functionality.

`carriercapture.core.transfer_coord.TransferCoordinate` ¶

Configuration coordinate for charge transfer calculations using Marcus theory.

Represents two diabatic potential energy surfaces (localized states) for calculating electron/hole transfer rates and charge carrier mobility.

This class implements Marcus theory as described in: Marcus, R. A. J. Chem. Phys. 24, 966 (1956)

Attributes:

Name	Type	Description
`name`	`str`	Identifier
`pot_1`	`Potential`	First diabatic state potential
`pot_2`	`Potential`	Second diabatic state potential
`Q_cross`	`float \| None`	Configuration coordinate at intersection point (amu^0.5·Å)
`E_cross`	`float \| None`	Energy at intersection point (eV)
`coupling`	`float \| None`	Electronic coupling Hab at intersection (eV)
`reorganization_energy`	`float \| None`	Marcus reorganization energy λ (eV)
`activation_energy`	`float \| None`	Classical activation barrier ΔG‡ (eV)
`transfer_rate`	`NDArray[float64] \| None`	Transfer rate vs temperature (s⁻¹)
`temperature`	`NDArray[float64] \| None`	Temperature grid (K)

Examples:

Basic usage for charge transfer:

>>> pot_1 = Potential.from_harmonic(hw=0.02, Q0=0.0, E0=0.0)
>>> pot_2 = Potential.from_harmonic(hw=0.02, Q0=5.0, E0=0.1)
>>> tc = TransferCoordinate(pot_1, pot_2, name="hole_transfer")
>>> tc.get_coupling()
>>> tc.get_reorganization_energy()
>>> tc.get_transfer_rate(temperature=np.linspace(100, 500, 50))

`init(pot_1, pot_2, name='')` ¶

Initialize a TransferCoordinate.

Parameters:

Name	Type	Description	Default
`pot_1`	`Potential`	First diabatic state potential	required
`pot_2`	`Potential`	Second diabatic state potential	required
`name`	`str`	Identifier for this transfer coordinate	`""`

`get_coupling(Q_cross=None)` ¶

Calculate electronic coupling between diabatic states.

The coupling Hab is half the energy splitting between adiabatic states at the intersection point of the diabatic surfaces.

Parameters:

Name	Type	Description	Default
`Q_cross`	`float`	Intersection point (amu^0.5·Å). If None, will find crossing automatically using find_crossing().	`None`

Returns:

Name	Type	Description
`coupling`	`float`	Electronic coupling Hab (eV)

Raises:

Type	Description
`ValueError`	If potentials don't have fit functions
`ValueError`	If no crossing point found

Notes

The adiabatic states are: E± = 0.5 * (E1 + E2) ± sqrt(0.25 * (E1 - E2)² + Hab²)

At the diabatic crossing (E1 = E2), the splitting is: E+ - E- = 2 * Hab

Examples:

>>> tc.get_coupling()  # Automatic crossing detection
0.015
>>> tc.get_coupling(Q_cross=2.5)  # Specify crossing point
0.018

`get_reorganization_energy()` ¶

Calculate Marcus reorganization energy.

The reorganization energy λ is the energy required to rearrange the nuclear coordinates from the equilibrium geometry of state 1 to that of state 2 (or vice versa), while remaining in state 1.

For symmetric potentials: λ = E1(Q2_min) - E1(Q1_min) = E2(Q1_min) - E2(Q2_min)

Returns:

Name	Type	Description
`lambda`	`float`	Reorganization energy (eV)

Raises:

Type	Description
`ValueError`	If potentials don't have fit functions

Notes

The reorganization energy is related to the Huang-Rhys factor S: λ = 2 * S * ℏω for harmonic potentials.

Examples:

>>> tc.get_reorganization_energy()
0.25

`get_activation_energy(delta_G=0.0)` ¶

Calculate classical activation energy for charge transfer.

Marcus formula for activation energy: ΔG‡ = (λ + ΔG)² / (4λ)

Where: λ = reorganization energy ΔG = reaction free energy (driving force)

For symmetric case (ΔG = 0): ΔG‡ = λ / 4

Parameters:

Name	Type	Description	Default
`delta_G`	`float`	Reaction free energy (eV) Positive = uphill, Negative = downhill	`0.0`

Returns:

Name	Type	Description
`barrier`	`float`	Activation energy (eV)

Raises:

Type	Description
`ValueError`	If reorganization energy hasn't been calculated

Examples:

>>> tc.get_activation_energy()  # Symmetric case
0.0625
>>> tc.get_activation_energy(delta_G=-0.1)  # Downhill
0.04

`get_transfer_rate(temperature, delta_G=0.0)` ¶

Calculate Marcus transfer rate vs temperature.

Marcus non-adiabatic transfer rate: k(T) = (2π/ℏ) * (1/sqrt(4πλkBT)) * Hab² * exp(-ΔG‡/kBT)

Where: ΔG‡ = (λ + ΔG)² / (4λ)

Parameters:

Name	Type	Description	Default
`temperature`	`NDArray[float64]`	Temperature grid (K)	required
`delta_G`	`float`	Reaction free energy (eV)	`0.0`

Returns:

Name	Type	Description
`rate`	`NDArray[float64]`	Transfer rate (s⁻¹)

Raises:

Type	Description
`ValueError`	If coupling or reorganization energy haven't been calculated

Examples:

>>> temps = np.linspace(100, 500, 50)
>>> tc.get_transfer_rate(temperature=temps)
array([1.23e12, 2.45e12, ...])

`calculate_mobility(temperature, distance, delta_G=0.0)` ¶

Calculate Einstein mobility from transfer rate.

Einstein relation: μ = e * D / (kB * T)

Where diffusion coefficient: D = d² * k / 2

For 1D hopping with distance d and rate k.

Parameters:

Name	Type	Description	Default
`temperature`	`NDArray[float64]`	Temperature grid (K)	required
`distance`	`float`	Hopping distance (Å)	required
`delta_G`	`float`	Reaction free energy (eV)	`0.0`

Returns:

Name	Type	Description
`mobility`	`NDArray[float64]`	Charge carrier mobility (cm²/(V·s))

Notes

This is a simplified model assuming 1D hopping between nearest-neighbor sites. Real materials require 3D treatment.

Examples:

>>> temps = np.linspace(100, 500, 50)
>>> tc.calculate_mobility(temperature=temps, distance=5.0)
array([1.2e-4, 2.3e-4, ...])

`to_dict()` ¶

Serialize transfer coordinate to dictionary.

Returns:

Name	Type	Description
`data`	`dict`	Dictionary representation

`from_dict(data)` `classmethod` ¶

Deserialize transfer coordinate from dictionary.

Parameters:

Name	Type	Description	Default
`data`	`dict`	Dictionary representation	required

Returns:

Name	Type	Description
`tc`	`TransferCoordinate`	Reconstructed transfer coordinate

Schrödinger Solver¶

Low-level functions for solving the 1D Schrödinger equation.

`carriercapture.core.schrodinger` ¶

1D Schrödinger equation solver using finite difference method and ARPACK.

This module provides functions to solve the time-independent Schrödinger equation for 1D potentials, which is critical for computing phonon states in carrier capture calculations.

Performance: This is the computational bottleneck (80-90% of runtime). Uses scipy's ARPACK wrapper (same FORTRAN backend as Julia) for efficient sparse eigenvalue solving.

`build_hamiltonian_1d(potential_func, Q)` ¶

Build sparse Hamiltonian matrix for 1D TISE using finite differences.

Uses three-point finite difference stencil for kinetic energy operator: T ≈ -ℏ²/(2m) * (ψ[i+1] - 2ψ[i] + ψ[i-1]) / h²

The Hamiltonian is: H = T + V where T is the kinetic energy (sparse tridiagonal) and V is the potential energy (diagonal).

Parameters:

Name	Type	Description	Default
`potential_func`	`callable`	Potential energy function V(Q) in eV Should accept numpy array and return numpy array	required
`Q`	`NDArray[float64]`	Configuration coordinate grid (amu^0.5·Å) Must be uniformly spaced	required

Returns:

Name	Type	Description
`H`	`csr_matrix`	Sparse Hamiltonian matrix (N × N) in CSR format

Notes

The finite difference method assumes a uniform grid spacing h = Q[1] - Q[0]. Boundary conditions are infinite potential walls (ψ = 0 at boundaries).

The mass is m = 1 amu, typical for nuclear configuration coordinates.

Examples:

>>> def harmonic(Q):
...     return 0.5 * 0.02 * Q**2
>>> Q = np.linspace(-10, 10, 1000)
>>> H = build_hamiltonian_1d(harmonic, Q)
>>> H.shape
(1000, 1000)

`normalize_wavefunctions(wavefunctions, grid_spacing)` ¶

Normalize wavefunctions so ∫|ψ|² dQ = 1.

Parameters:

Name	Type	Description	Default
`wavefunctions`	`NDArray[float64]`	Wavefunctions, shape (nev, N_Q) Each row is a wavefunction	required
`grid_spacing`	`float`	Spacing between grid points (amu^0.5·Å)	required

Returns:

Name	Type	Description
`normalized`	`NDArray[float64]`	Normalized wavefunctions, same shape as input

Notes

Uses trapezoidal rule for numerical integration: ∫|ψ|² dQ ≈ h * Σ|ψ|²

Examples:

>>> wf = np.random.rand(10, 1000)  # 10 wavefunctions, 1000 points
>>> h = 0.01
>>> wf_norm = normalize_wavefunctions(wf, h)
>>> np.allclose(np.sum(wf_norm**2, axis=1) * h, 1.0)
True

`solve_schrodinger_1d(potential_func, Q, nev=30, maxiter=None)` ¶

Solve 1D time-independent Schrödinger equation using ARPACK.

Uses scipy.sparse.linalg.eigsh (which wraps ARPACK FORTRAN library) to find the lowest nev eigenvalues and eigenvectors of the Hamiltonian.

This is the computational bottleneck (80-90% of total runtime) for carrier capture calculations.

Parameters:

Name	Type	Description	Default
`potential_func`	`callable`	Potential energy function V(Q) in eV Should accept numpy array and return numpy array	required
`Q`	`NDArray[float64]`	Configuration coordinate grid (amu^0.5·Å) Must be uniformly spaced, typically 500-5000 points	required
`nev`	`int`	Number of eigenvalues to compute (lowest energy states) Typical range: 30-180	`30`
`maxiter`	`int`	Maximum number of ARPACK iterations If None, defaults to max(nev * len(Q), 1000)	`None`

Returns:

Name	Type	Description
`eigenvalues`	`NDArray[float64]`	Phonon energies (eV), shape (nev,) Sorted in ascending order
`eigenvectors`	`NDArray[float64]`	Normalized wavefunctions, shape (nev, len(Q)) Each row is a wavefunction ψ_n(Q)

Raises:

Type	Description
`RuntimeError`	If ARPACK fails to converge within maxiter iterations
`ValueError`	If Q grid is not uniformly spaced

Notes

Grid must be uniform for finite difference method
Wavefunctions are normalized: ∫|ψ|² dQ = 1
Uses mass m = 1 amu (typical for nuclear coordinates)
Boundary conditions: ψ = 0 at Q_min and Q_max (infinite walls)

Performance

N=5000, nev=60: ~0.5-1s (comparable to Julia)
Bottleneck is ARPACK eigenvalue solver
Same FORTRAN backend as Julia, so performance is identical

Examples:

Harmonic oscillator:

>>> def harmonic(Q):
...     hw = 0.02  # ℏω = 20 meV
...     return 0.5 * hw * Q**2
>>> Q = np.linspace(-20, 20, 5000)
>>> eigenvalues, eigenvectors = solve_schrodinger_1d(harmonic, Q, nev=10)
>>> eigenvalues[:3]  # First 3 eigenvalues
array([0.01, 0.03, 0.05])  # E_n = ℏω(n + 1/2)

Constants¶

Physical constants used throughout CarrierCapture.

`carriercapture._constants` ¶

Physical constants used in carrier capture calculations.

All constants are in SI-compatible units suitable for semiconductor physics: - Energy: eV (electronvolts) - Length: Å (angstroms) or cm - Mass: amu (atomic mass units) - Temperature: K (kelvin)

`OCC_CUTOFF = 1e-05` `module-attribute` ¶

Maximum occupation of highest eigenvalue for convergence.

If occ(ε_max) >= OCC_CUTOFF, the partition function is not converged and more eigenvalues (nev) should be computed.

Core API¶

Potential¶

carriercapture.core.potential.Potential ¶

__init__(name='', Q0=0.0, E0=0.0, nev=30, temperature=293.15, Q_data=None, E_data=None) ¶

from_file(filepath, name='', resolution=3001, **kwargs) classmethod ¶

from_harmonic(hw, Q0=0.0, E0=0.0, Q_range=(-20.0, 20.0), npoints=5000, **kwargs) classmethod ¶

fit(fit_type='spline', **fit_kwargs) ¶

filter_thermally_accessible(thermal_energy=None) ¶

solve(nev=None, maxiter=None) ¶

__call__(Q) ¶

copy() ¶

to_dict() ¶

from_dict(data) classmethod ¶

ConfigCoordinate¶

carriercapture.core.config_coord.ConfigCoordinate ¶

__init__(pot_i, pot_f, name='', W=0.0, degeneracy=1) ¶

calculate_overlap(Q0, cutoff=0.25, sigma=0.025) ¶

calculate_capture_coefficient(volume, temperature) ¶

to_dict() ¶

from_dict(data) classmethod ¶

TransferCoordinate¶

carriercapture.core.transfer_coord.TransferCoordinate ¶

__init__(pot_1, pot_2, name='') ¶

get_coupling(Q_cross=None) ¶

get_reorganization_energy() ¶

get_activation_energy(delta_G=0.0) ¶

get_transfer_rate(temperature, delta_G=0.0) ¶

calculate_mobility(temperature, distance, delta_G=0.0) ¶

to_dict() ¶

from_dict(data) classmethod ¶

Schrödinger Solver¶

carriercapture.core.schrodinger ¶

build_hamiltonian_1d(potential_func, Q) ¶

normalize_wavefunctions(wavefunctions, grid_spacing) ¶

solve_schrodinger_1d(potential_func, Q, nev=30, maxiter=None) ¶

Constants¶

carriercapture._constants ¶

OCC_CUTOFF = 1e-05 module-attribute ¶

`carriercapture.core.potential.Potential` ¶

`init(name='', Q0=0.0, E0=0.0, nev=30, temperature=293.15, Q_data=None, E_data=None)` ¶

`from_file(filepath, name='', resolution=3001, **kwargs)` `classmethod` ¶

`from_harmonic(hw, Q0=0.0, E0=0.0, Q_range=(-20.0, 20.0), npoints=5000, **kwargs)` `classmethod` ¶

`fit(fit_type='spline', **fit_kwargs)` ¶

`filter_thermally_accessible(thermal_energy=None)` ¶

`solve(nev=None, maxiter=None)` ¶

`call(Q)` ¶

`copy()` ¶

`to_dict()` ¶

`from_dict(data)` `classmethod` ¶

`carriercapture.core.config_coord.ConfigCoordinate` ¶

`init(pot_i, pot_f, name='', W=0.0, degeneracy=1)` ¶

`calculate_overlap(Q0, cutoff=0.25, sigma=0.025)` ¶

`calculate_capture_coefficient(volume, temperature)` ¶

`to_dict()` ¶

`from_dict(data)` `classmethod` ¶

`carriercapture.core.transfer_coord.TransferCoordinate` ¶

`init(pot_1, pot_2, name='')` ¶

`get_coupling(Q_cross=None)` ¶

`get_reorganization_energy()` ¶

`get_activation_energy(delta_G=0.0)` ¶

`get_transfer_rate(temperature, delta_G=0.0)` ¶

`calculate_mobility(temperature, distance, delta_G=0.0)` ¶

`to_dict()` ¶

`from_dict(data)` `classmethod` ¶

`carriercapture.core.schrodinger` ¶

`build_hamiltonian_1d(potential_func, Q)` ¶

`normalize_wavefunctions(wavefunctions, grid_spacing)` ¶

`solve_schrodinger_1d(potential_func, Q, nev=30, maxiter=None)` ¶

`carriercapture._constants` ¶

`OCC_CUTOFF = 1e-05` `module-attribute` ¶