A modular and extensible Python package for calculating electronic properties of graphene and multilayer graphene systems using tight-binding theory.

This package provides a clean, reusable framework for simulating graphene properties using the tight-binding approximation. It's organized to facilitate easy addition of new properties and calculations in the future.

• Calculate band structure around K-point
• Determine minimum direct band gap
• Support for monolayer to bulk graphite (configurable layers)
• Multiple stacking patterns: AB, ABC, ABA

• Complex dielectric function: ε(ω) = ε₁(ω) + iε₂(ω)
• Refractive index and extinction coefficient
• Reflectance and transmittance at specific wavelengths
• Optical conductivity via Kubo formula

• Electron energy loss spectroscopy (EELS) simulations
• Drude-Lorentz oscillator model for realistic lineshapes
• π and π+σ plasmon peaks
• Layer-dependent oscillator parameters

• Broadband optical response across 0.5–6 eV
• Multi-panel visualization of n, k, R, T, and absorbance
• Layer-vs-energy heatmaps for absorbance, reflectance, and transmittance

• Displacement-field control of optical constants
• Quick comparison utilities for neutral vs biased stacks
• Summary metrics (Δn, ΔA) at target wavelengths

• Gaussian-broadened DOS from tight-binding spectra
• Layer-resolved evolution from graphene to graphite
• Dirac cone sampling utilities for linear dispersion visualization

graphene_simulation/
├── graphene_tb/
│   ├── __init__.py           # Package entry point
│   ├── config.py             # Configuration classes & constants
│   ├── hamiltonian.py        # TB Hamiltonian & band structure
│   ├── optical.py            # Dielectric function & optical properties
│   ├── eels.py              # EELS calculations
│   └── plotting.py          # Visualization utilities
├── Demo_Graphene_TB_Simulation.ipynb  # Comprehensive demonstration
└── README.md                 # This file

Defines physical constants and the GrapheneConfig dataclass for system configuration.

Key Classes:

• GrapheneConfig: Configurable system parameters (hopping integrals, layer numbers, applied fields)

Physical Constants:

• C-C bond length: 0.142 nm
• Interlayer distance: 0.335 nm
• K-point coordinates

Core tight-binding calculations for band structure.

Key Functions:

• build_full_hamiltonian(): Construct multilayer Hamiltonian
• find_band_gap(): Calculate minimum band gap
• calculate_band_structure(): Full band structure calculation
• get_stacking_sequence(): Handle different stacking patterns

Optical response calculations.

Key Functions:

• calculate_dielectric_function(): ε(ω) using Kubo formula + Kramers-Kronig
• calculate_refractive_index(): Convert ε to n and k
• calculate_interband_optical_conductivity(): σ(ω) from band transitions
• compute_optical_response(): Full spectral response (n, k, R, T, A, σ)
• compute_layer_optical_map(): Layer-vs-energy response grids

Electron energy loss spectroscopy.

Key Functions:

• epsilon_drude_lorentz(): Drude-Lorentz dielectric model
• graphene_oscillators_for_layers(): Layer-dependent oscillator parameters
• loss_function_from_eps(): EELS loss function = Im[-1/ε]

Electronic density of states and Dirac cone helpers.

Key Functions:

• calculate_density_of_states(): Gaussian-broadened DOS over energy grid
• sample_dirac_cone(): High-resolution sampling around the K-point for Dirac cone analysis

Publication-quality visualization utilities.

Key Functions:

• plot_band_gap_vs_layers(): Band gap evolution
• plot_band_structure_path(): High-symmetry band structure lines
• plot_optical_constants(): n and k with dual y-axes
• plot_reflectance_transmittance(): R and T with dual y-axes
• plot_optical_response_panels(): Multi-panel spectral overview (n, k, R, T, A, α)
• plot_optical_phase_map(): Heatmaps for layer-energy response
• plot_density_of_states(): Overlay DOS curves for multiple systems
• plot_dirac_cone(): 3D and contour visualization of Dirac dispersion
• set_academic_style(): Apply academic formatting

import sys
from pathlib import Path
sys.path.insert(0, str(Path('.').resolve()))

import graphene_tb
from graphene_tb import GrapheneConfig, find_band_gap

# Create configuration
config = GrapheneConfig(N_layers=2, stacking="AB")

# Calculate band gap
gap, k_point = find_band_gap(config)
print(f"Band gap: {gap:.4f} eV at k = {k_point}")

See Demo_Graphene_TB_Simulation.ipynb for a comprehensive example showing:

1. Band gap calculations vs layer number
2. Optical properties (n, k, R, T)
3. EELS spectra with Drude-Lorentz model
4. Comparison between monolayer, few-layer, and graphite

from graphene_tb import compute_optical_response, compute_layer_optical_map
from graphene_tb.plotting import plot_optical_response_panels, plot_optical_phase_map

energies = np.linspace(0.5, 6.0, 240)
config = GrapheneConfig(N_layers=4)

# Full spectral response (n, k, R, T, A, α)
response = compute_optical_response(config, energies, k_sampling=60, broadening=0.08)
plot_optical_response_panels(energies, [response], ["Few-layer (N=4)"])

# Layer-energy absorbance map for N = 1…20
layers = np.arange(1, 21)
absorbance_map, responses = compute_layer_optical_map(
    config,
    energies,
    layers,
    quantity="absorbance",
    k_sampling=45,
)
plot_optical_phase_map(energies, layers, absorbance_map, "Absorbance")
plt.show()

from graphene_tb import calculate_density_of_states, sample_dirac_cone
from graphene_tb.plotting import plot_density_of_states, plot_dirac_cone

energy_grid = np.linspace(-3.0, 3.0, 400)
configs = [
    ("Graphene", GrapheneConfig(N_layers=1)),
    ("Few-layer (N=4)", GrapheneConfig(N_layers=4)),
    ("Graphite-like", GrapheneConfig(N_layers=20)),
]

dos_curves = [
    calculate_density_of_states(cfg, energy_grid, k_sampling=70, broadening=0.04)
    for _, cfg in configs
]
plot_density_of_states(energy_grid, dos_curves, [label for label, _ in configs])

kx_mesh, ky_mesh, valence, conduction = sample_dirac_cone(configs[0][1], resolution=61)
plot_dirac_cone(kx_mesh, ky_mesh, valence, conduction)
plt.show()

• γ₀: In-plane nearest-neighbor hopping = 3.0 eV
• γ₁: Interlayer hopping (A-B) = 0.39 eV
• γ₃, γ₄: Higher-order corrections (currently unused)

• AB: Alternating ABA... stacking (most common)
• ABC: Rhombohedral stacking
• ABA: Mirror symmetry

• field_strength: Perpendicular displacement field (V/nm)
• screening_factor: Field screening (0=none, 1=full)

• Monolayer: 0 eV (semi-metallic, Dirac point at K)
• Bilayer: ~0 eV (slight opening with symmetry breaking)
• Few-layer: Rapidly vanishing
• Graphite (N→∞): ~0 eV (metallic, overlapping bands)

• Monolayer: n ≈ 1.63, highly transparent
• 4-layer: n ≈ 5-6, semi-transparent
• Graphite: n ≈ 15-17, opaque metallic
• Layer-energy heatmap: peak absorbance concentrates near 4.5–5.5 eV for 6–12 layers
• Electric field (0.6 V/nm) on N=4 stack boosts n and absorbance while suppressing transmittance

• Graphene retains linear DOS near the Dirac point, vanishing at E=0
• Intermediate layers broaden DOS and introduce shoulder features
• Graphite limit shows finite DOS at E=0, signalling semimetallic overlap

• π plasmon: 4.7 eV (monolayer) → 6.3 eV (graphite)
• π+σ plasmon: ~14.8-15.4 eV (weak layer dependence)

1. Create new module (e.g., density_of_states.py):

# graphene_tb/density_of_states.py
from .hamiltonian import build_full_hamiltonian
from .config import GrapheneConfig

def calculate_dos(config, energy_range, n_bins=200):
    # Implementation here
    pass

2. Add to __init__.py:

from .density_of_states import calculate_dos
__all__ = [..., 'calculate_dos']

3. Add plotting function to plotting.py

4. Use in notebook:

from graphene_tb import calculate_dos
dos = calculate_dos(config, energy_range)

• Density of States (DOS): Via tetrahedron method or histogram
• Electronic Transport: Conductivity, mobility, Hall effect
• Phonon Properties: Dispersion relations, Raman spectroscopy
• Magnetic Properties: With applied magnetic field
• Berry Phase & Topology: Topological properties
• STM/STS: Tunneling spectroscopy simulations

• Band gap calculations: O(N²) in k-mesh size
• Optical properties: O(N²·M) where M = energy points
• EELS: Similar to optical properties
• Typical runtime: <1 min for full parameter sweep on standard laptop

• Reduce k_mesh_size for quick scans
• Use k_sampling parameter in optical calculations
• Vectorize energy loops where possible

This implementation uses simplified tight-binding with:

• ✓ Correct band structure topology
• ✓ Realistic plasmon energies
• ✓ Qualitative layer dependence
• ✗ Doesn't include excitonic effects
• ✗ No many-body corrections
• ✗ Limited to K-point vicinity

For publication-quality results, compare with:

• First-principles calculations (DFT)
• Experimental EELS data
• Published graphene optical measurements

• Tight-Binding Theory: Slonczewski & Weiss, Phys. Rev. 109, 272 (1958)
• Graphene Band Structure: Wallace, Phys. Rev. 71, 622 (1947)
• Optical Properties: Nair et al., Science 320, 1308 (2008)
• EELS Spectroscopy: Egerton, Reports on Progress in Physics 72, 016502 (2009)

This package is provided for educational and research purposes.

Graphene Simulation Team

v1.2.0 - Added DOS & Dirac cone utilities with dedicated plotting and documentation v1.1.0 - Added broadband optical response utilities, layer-energy maps, and field-tunable optics visuals v1.0.0 - Initial release with band structure, optical properties, and EELS

Last Updated: October 2025