almaBTE: a solver of the space-time dependent Boltzmann transport equation for phonons in structured materials

Jesús Carrete; Bjorn Vermeersch; Ankita Katre; Ambroise van Roekeghem; Tao Wang; Georg K. H. Madsen; Natalio Mingo

Source-linked AI summary

almaBTE: a solver of the space-time dependent Boltzmann transport equation for phonons in structured materials

Jesús Carrete, Bjorn Vermeersch, Ankita Katre, Ambroise van Roekeghem, Tao Wang, Georg K. H. Madsen, Natalio Mingo

arXiv:1704.04142v1Apr 13, 2017physics.comp-ph

TL;DR

Thermal transport in technologically relevant semiconductor structures requires methods beyond bulk models because phonon flow can be multiscale and spatially complex. almaBTE provides first-principles, space- and time-dependent BTE solvers and supporting models for bulk materials, alloys, thin films, and superlattices. The package offers a modular and flexible framework for predictive thermal-transport analysis from the nano- to the microscale, while specific approximations constrain some regimes.

Problem
Structured semiconductor heat transport can involve simultaneous ballistic and diffusive phonon flow, while position-dependent BTE calculations face a six-dimensional curse of dimensionality.
Method
almaBTE combines first-principles phonon data with analytical, Monte Carlo, bulk, thin-film, alloy, and superlattice BTE solvers in a modular software package.
Results
almaBTE provides a richer and more accurate picture of thermal transport from the nano- to the microscale than classical equations and continuum models.
Takeaways & Limitations
The open-source package supports predictive studies of phonon-mediated heat transport across bulk and structured semiconducting systems.
Takeaways & Limitations
The cross-plane suppression function is not rigorously exact, although prior work found it reproduces effective-conductivity Monte Carlo simulations within a few percent.

Abstract

from arXiv · show

almaBTE is a software package that solves the space- and time-dependent Boltzmann transport equation for phonons, using only ab-initio calculated quantities as inputs. The program can predictively tackle phonon transport in bulk crystals and alloys, thin films, superlattices, and multiscale structures with size features in the nm-$μ$m range. Among many other quantities, the program can output thermal conductances and effective thermal conductivities, space-resolved average temperature profiles, and heat-current distributions resolved in frequency and space. Its first-principles character makes almaBTE especially well suited to investigate novel materials and structures. This article gives an overview of the program structure and presents illustrative examples for some of its uses.

PROGRAM SUMMARY

almaBTE calculates temperature profiles, thermal flux distributions, and effective thermal conductivities in structured systems where phonons carry heat.

almaBTE is a C++ software package for calculating temperature profiles, thermal flux distributions, and effective thermal conductivities in phonon-mediated structured systems.It runs on Linux and macOS and uses ab-initio-derived inputs with supporting libraries including BOOST, MPI, Eigen, HDF5, and spglib.

1. Introduction

Thermal management in semiconductor technologies requires modeling phonon transport across multiscale, structured systems. almaBTE addresses the resulting space- and time-dependent BTE challenge with a first-principles framework combining multiple solvers and models.

Phonons are the main heat carriers in semiconductors and insulators, making multiscale vibration physics central to heat dissipation.
A 50 K reduction in an LED or transistor’s active-region temperature can increase time to failure by one order of magnitude.
Position-dependent phonon BTE calculations face a curse of dimensionality because they involve six spatial dimensions.These dimensions combine three real-space and three reciprocal-space coordinates, making direct linear-system solutions unfeasible.
almaBTE provides space- and time-dependent phonon BTE solvers and auxiliary programs for first-principles studies of structured semiconducting systems.
The package includes improved bulk functionality, thin-film conductivity models, a variance-reduced one-dimensional Monte Carlo solver, superlattice models, and an analytical one-dimensional solver.

2. General structure of almaBTE

almaBTE uses a modular workflow in which effective-crystal data are generated from ab-initio inputs, stored in HDF5, and passed to property explorers and BTE solvers. Its C++ implementation supports portable, parallel execution.

A from-scratch calculation begins by creating effective crystals that serve as building blocks for complex structures.Each effective crystal can be assigned an atomic structure and phonon spectrum.
Material descriptions are stored in HDF5 groups for crystal structure, reciprocal-space phonon data, elastic scattering, and three-phonon processes.The three-phonon group contains the processes and amplitudes needed to compute intrinsic scattering rates.
HDF5 files can represent single crystals, virtual-crystal-approximation alloys, and superlattices, using ab-initio geometries and second- and third-order potential-energy derivatives.
The resulting HDF5 files feed thermal-property explorers and BTE solvers, each configured through its own XML input file.Outputs use text-based formats readable and parseable by common data-analysis software.
almaBTE is written in standard C++11, supports Linux and macOS, and uses MPI parallelization for performance-critical calculations.The reported performance scales almost linearly on non-I/O-constrained calculations over at least the tested processor range.

3. Programs included in the package

The package combines builders for crystals, alloys, and superlattices with analytical and numerical transport solvers. Its models use phonon spectra, scattering processes, composition profiles, and structured perturbations derived from first-principles inputs.

3.1. VCAbuilder: VCAbuilder creates HDF5 descriptions of single crystals and virtual-crystal-approximation alloys.
3.1. VCAbuilder: For single crystals, VCAbuilder constructs reciprocal-space dynamical matrices from structural data and harmonic interatomic force constants.It can also incorporate information for polar bonds.
3.1. VCAbuilder: VCAbuilder searches allowed energy- and momentum-conserving three-phonon processes and computes their scattering amplitudes from anharmonic interatomic force constants.
3.1. VCAbuilder: Symmetry reductions replace equivalent wavevectors with irreducible wavevectors to manage the 10^6–10^12 allowed processes typical of reciprocal-space grids.
3.1. VCAbuilder: For alloys, effective parameters are composition-weighted averages of component properties, avoiding additional ab-initio calculations for new concentrations.The approximation does not account for correlations, local relaxations, or alloying-induced electronic-structure changes.
3.2. superlattice builder: The superlattice builder models alternating compounds with layer-dependent compositions, including digital and segregated profiles and their one-dimensional versus three-dimensional disorder.
3.2. superlattice builder: Superlattice scattering combines Tamura mass-disorder scattering with barrier scattering from the average mass profile along the growth direction.The barrier treatment uses a causal Green’s function and t matrix while neglecting force-constant changes within the virtual-crystal approximation.
3.2. superlattice builder: Cubic band interpolation with continuous derivatives removes artificial divergences introduced where piecewise-linear band segments meet with different slopes.The resulting analytical sub-segment contributions reduce the number of grid points needed for the one-dimensional Green’s function.

3.3. kappa Tsweep

almaBTE efficiently evaluates bulk thermal conductivity across ambient temperatures by reusing stored phonon-scattering information and computing temperature-dependent rates as needed.

The temperature-sweep executable computes thermal conductivity κ(T) efficiently across ambient temperatures.Stored three-phonon processes and scattering matrix elements do not need recomputation.
Temperature-independent two-phonon scattering rates are precomputed, while three-phonon scattering rates are evaluated separately at each temperature.
Mass-disorder scattering contributes to two-phonon rates for crystals containing multiple isotopes or alloying.
Additional elastic-scattering datasets, such as dislocations, can be combined with mass-disorder scattering, whereas superlattices replace this treatment with layer-resolved rates.
Under the relaxation time approximation, thermal conductivity is obtained from mode heat capacity, group velocity, and mean free path contributions.The mean free path is Λ(T) = |v|τ(T), with sums over phonon branches and Brillouin-zone averages.
Beyond RTA, almaBTE solves the linearized BTE directly after reducing unknowns to irreducible wavevector points using crystal rotational symmetries.Eigen routines are used for the resulting system.

3.4. kappa crossplanefilms, kappa inplanefilms

almaBTE models thin-film thermal transport in both in-plane and cross-plane configurations using suppression functions that account for boundary scattering and phonon direction.

The thin-film executables assess one-dimensional thermal transport for both in-plane and cross-plane configurations.
In-plane transport: The in-plane suppression function is evaluated mode by mode in wavevector space, preserving crystal anisotropy rather than assuming isotropic dispersions.
In-plane transport: The suppression function accounts for additional boundary scattering, while the angle ϑ measures group velocity relative to the transport axis.
In-plane transport: For in-plane transport, the geometry uses the film normal n and an orthogonal transport vector u.
In-plane transport: Film boundaries share a frequency-independent specularity p, and coordinates are rotated so the film normal aligns with the transformed z axis.
In-plane transport: In the thick-film limit, S∥(K →0) = 1 regardless of p, recovering bulk transport; modes parallel to the film contribute their full nominal conductivity.
Cross-plane transport: Cross-plane transport is specified by u ≡ n and assumes perfectly absorbing black-body film boundaries.
Cross-plane transport: For cross-plane transport, K = Λ|cos ϑ|/L avoids coordinate transforms, while the approximate suppression function reproduces Monte Carlo conductivities within a few percent.

3.5. cumulativecurves

almaBTE produces cumulative curves showing how phonon-mode contributions to bulk heat capacity and thermal conductivity accumulate when resolved by selected transport parameters.

The cumulative-curve executable reveals phonon-mode contributions to bulk heat capacity and thermal conductivity.
Curves can be resolved by mean free path, relaxation time, projected mean free path, frequency, angular frequency, or energy.The projected mean free path Λ|cos θ| is central to cross-plane film transport.

3.6. transient analytic1d

The analytic1d module provides semi-analytic solutions for one-dimensional, time-dependent phonon transport under the relaxation time approximation. It computes temperature profiles and related observables by semi-analytically or numerically inverting transformed-domain solutions.

3.6. transient analytic1d: analytic1d solves the single-pulse response of the RTA-BTE for one-dimensional transient phonon transport in infinite bulk media.The thermal field depends only on the Cartesian coordinate x.
3.6.1. Temperature profiles ∆T(x, t): A planar source at x = 0 injects 1 J/m2 at t = 0, with energy distributed among phonon modes according to heat-capacity contributions.
3.6.1. Temperature profiles ∆T(x, t): The module obtains real-space and time quantities by inverting Fourier-Laplace solutions, using piecewise second-order Taylor expansions for Fourier inversion.The expanded intervals are integrated analytically.
3.6.1. Temperature profiles ∆T(x, t): The simplified inversion assumes weakly quasiballistic conditions and ignores purely ballistic effects, while performing well beyond characteristic relaxation times of typically 1–10 ns.
3.6.1. Temperature profiles ∆T(x, t): The desired temperature profile is obtained by dividing the deviational energy density by the bulk heat capacity.The source temperature rise can also be evaluated at x = 0 using linear quadrature.
3.6.3. Mean square displacement σ2(t): The thermal mean-square displacement is evaluated through Laplace-domain moments and numerically inverted with the Gaver-Stehfest method.The resulting expression applies across ballistic, quasiballistic, and diffusive time scales.

3.7. steady montecarlo1d

The steady montecarlo1d module simulates one-dimensional steady transport in multilayer structures using first-principles, wavevector-resolved phonon properties. Variance-reduced deviational particles are launched, propagated through scattering and interfaces, and used to estimate temperature profiles, heat flux, and effective thermal metrics.

3.7. steady montecarlo1d: steady montecarlo1d directly simulates steady-state thermal transport in one-dimensional multilayered structures with isothermal boundaries.
3.7. steady montecarlo1d: The method uses wavevector-resolved, first-principles phonon properties rather than parameterized dispersions or frequency-only formulations, supporting anisotropic crystals and predictive modeling of novel materials.
3.7. steady montecarlo1d: Deviational particles represent fixed amounts of power in the steady-state formulation, exploiting time-translation invariance to avoid an equilibration phase.
3.7. steady montecarlo1d: Under the linearized regime, particles propagate independently with phonon properties fixed at Tref = (Thot + Tcold)/2, enabling parallelization without spatial or temporal discretization.
3.7. steady montecarlo1d: Particles undergo intrinsic scattering, interface transmission or reflection, and reservoir absorption as trajectories encounter material boundaries and isothermal reservoirs.
3.7.6. Evaluation of temperature profile: Temperature profiles are accumulated from trajectory subsegments crossing spatial bins, whereas heat flux is obtained by accumulating signed particle contributions at bin centers.
3.7.7. Effective thermal metrics: The simulated net heat flux yields apparent thermal conductivity and effective resistivity or conductance for the complete multilayer structure.

4. Examples of application

The examples apply almaBTE to bulk, thin-film, and layered systems, combining first-principles calculations with analytical and Monte Carlo transport methods. They demonstrate quantitative conductivity predictions, treatment of complex β-Ga2O3, and non-Fourier effects in Si/Ge bilayers.

4.1. Thermal conductivity of bulk materials: RTA conductivities agree with full BTE results within ∼5% for most common semiconductors, except diamond where normal scattering requires proper treatment.The comparison includes diamond, Si, and wurtzite GaN.
4.2. Effective thermal conductivity of thin films: β-Ga2O3 provides a challenging thin-film test case because its monoclinic structure is complex and strongly thermally anisotropic.The conventional crystal cell contains 20 atoms, while the primitive cell contains 10 atoms.
4.2. Effective thermal conductivity of thin films: The β-Ga2O3 thin-film calculation uses a 12 × 12 × 12 wavevector grid and produces an HDF5 file exceeding 8 GiB because of numerous three-phonon processes.The low crystal symmetry drives the large number of processes considered.
4.2. Effective thermal conductivity of thin films: 21 W/m/K is the calculated room-temperature bulk conductivity along β-Ga2O3’s most conductive (010) axis, agreeing with experimental measurements.The result differs from the 27 W/m/K value reported elsewhere and is consistent with prior ShengBTE calculations except along (001).
4.3. Cumulative conductivity: Cumulative conductivity curves resolve mean free path and energy contributions for wurtzite GaN, AlN, and Al0.5Ga0.5N.The calculations use 24 × 24 × 24 wavevector grids.
4.5. Steady-state Si/Ge bilayer in the steady state: In a 200 nm Si/Ge bilayer, temperature profiles become non-linear and show reservoir discontinuities from ballistic contact resistance, alongside a substantial Si/Ge interfacial drop.The spectral heat-flux map attributes the interface resistance to severe phonon-frequency mismatch: optical modes reflected from Si must scatter to lower frequencies before entering Ge.

5. Summary and conclusions

The paper positions first-principles thermal transport as predictive for crystalline semiconductors but insufficient alone for technologically relevant structures with defects and interfaces. almaBTE extends this capability across nanoscale-to-microscale systems through modular analytical and Monte Carlo BTE solvers.

5. Summary and conclusions: First-principles thermal-transport methods provide predictive thermal-conductivity estimates for crystalline semiconductors.The paper notes that these methods have recently emerged as powerful tools.
5. Summary and conclusions: Single-crystal conductivity is an upper bound for actual samples because crystallographic defects and interfaces affect phonon scattering.This scope boundary motivates treatment of more realistic structured systems.
5. Summary and conclusions: almaBTE combines modular analytical and Monte Carlo phonon-BTE approaches to model heat transport from the nano- to the microscale across a wider range of systems.The package is open source and accompanied by a database of publication-quality input files.

A. Cubic interpolation approach to obtaining the 1D Green’s function

The appendix develops a cubic-interpolation method for evaluating the 1D Green’s function by sampling dispersion-related quantities on a regular Brillouin-zone partition. It converts the integral into weighted sums and evaluates the resulting rational integrals through separate treatments of the denominator’s root structures.

A. Cubic interpolation approach to obtaining the 1D Green’s function: The 1D Green’s-function contribution from each phonon branch is expressed as an angular-frequency-dependent integral.The branch dispersion is represented by E(ϕ) = ω2(ϕ), with m(ϕ) as a smooth function.
A. Cubic interpolation approach to obtaining the 1D Green’s function: A regular partition of the 1D Brillouin zone samples m, E, and E′ at grid points before interpolation within each subinterval.The grid uses ϕn = 2πn/N.
A. Cubic interpolation approach to obtaining the 1D Green’s function: Cubic interpolants recast the Green’s-function integral as a weighted sum of sampled values of m(ϕ), with cyclic indexing across the grid.The partial weights are obtained from rational integrals.
A. Cubic interpolation approach to obtaining the 1D Green’s function: The rational-integral evaluation treats separately the cases of three real denominator roots and one real root.Reduced coefficients and discriminants determine the relevant root structure and intermediate quantities used to construct the weights.
A. Cubic interpolation approach to obtaining the 1D Green’s function: The method computes real and imaginary weight components using indicator and auxiliary functions after determining the appropriate roots and intermediate quantities.The displayed development defines functions including F1 and F3 for the weight formulas.