This notebook can be run on Google Colab.
In Colab, you can enable the GPU acceleration from Edit
> Notebook Settings
> Accelerator
> GPU
.
Amorphous silicon tutorial¶
pip install
ALDo¶
[ ]:
%%capture
! pip install kaldo
Fetch supplyment data remotely¶
[ ]:
%%capture
# Activate large file transfer in github
! git lfs install --skip-repo --skip-smudge
# Git clone kALDo
! git clone https://github.com/nanotheorygroup/kaldo.git
# Copy si-amorphous structure and force constant files
! cp -r /content/kaldo/kaldo/tests/si-amorphous ./
! rm -rf kaldo sample_data
Thermal transport simulation for Amorphous-Silicon (a-Si)¶
[ ]:
from kaldo.forceconstants import ForceConstants
# Read amporphous silicon with eskm format
forceconstants = ForceConstants.from_folder(folder='si-amorphous', format='eskm')
Create phonons object¶
[ ]:
from kaldo.phonons import Phonons
# Create a phonon object
# is_classic flag can be true or false, as
# well as 0 (quantum) or 1 (classic)
phonons = Phonons (forceconstants=forceconstants,
is_classic=0,
temperature=300,
folder='si-amorphous',
third_bandwidth=0.5/4.136,
broadening_shape='triangle',
storage='numpy')
Access and visualize properties calculated during simulations¶
[ ]:
import matplotlib.pyplot as plt
plt.style.use('seaborn-v0_8')
# Direct access to properties
# calculated during the simulation
frequency = phonons.frequency.flatten(order='C')
bandwidth = phonons.bandwidth.flatten(order='C')
# Plot phonon bandwidth vs frequency
print('\n')
plt.plot(frequency[3:],bandwidth[3:],'.',markersize=10,label= 'broadening shape: ' + str(phonons.broadening_shape))
plt.ylabel('$\Gamma$ (THz)', fontsize=25, fontweight='bold')
plt.xlabel("$\\nu$ (THz)", fontsize=25, fontweight='bold')
plt.ylim([bandwidth.min(), bandwidth.max()])
plt.legend(loc=2,frameon = False)
plt.grid()
plt.show()
phase_space = phonons.phase_space.flatten(order='C')
# Plot phase space vs frequency
print('\n')
plt.figure()
plt.plot(frequency[3:], phase_space[3:], '.', markersize=10)
plt.ylabel ("Phase space", fontsize=25, fontweight='bold')
plt.xlabel("$\\nu$ (THz)", fontsize=25, fontweight='bold')
plt.grid()
plt.show()
Calculate and visualize
and
¶
[ ]:
from kaldo.conductivity import Conductivity
from kaldo.controllers import plotter
import numpy as np
# Compute conductivity with per phonon mode basis using qhgk method
kappa_qhgk_per_mode = np.einsum('maa->m',1/3*Conductivity(phonons=phonons, method='qhgk',n_iterations=20).conductivity)
# Compute cumulative conductivity by frequency with qhgk method
kappa_qhgk_cum_freq = plotter.cumulative_cond_cal(phonons.frequency,Conductivity(phonons=phonons, method='qhgk').conductivity,phonons.n_phonons)
print('\n')
# Visualize the per mode conductivity vs frequency
plt.figure()
plt.plot(frequency[3:],kappa_qhgk_per_mode[3:],'.',label='QHGK',ms=8)
plt.xlabel ("$\\nu$ (THz)", fontsize=25, fontweight='bold')
plt.ylabel(r'$\kappa(W/m/K)$', fontsize=25, fontweight='bold')
plt.legend(loc=1,frameon=False)
plt.grid()
plt.show()
print('\n')
# Visualize the cumulative conductivity vs frequency
plt.figure()
plt.plot(frequency[3:],kappa_qhgk_cum_freq[3:],'.',label='QHGK',ms=8)
plt.xlabel ("$\\nu$ (THz)", fontsize=25, fontweight='bold')
plt.ylabel(r'$\kappa_{cum}(W/m/K)$', fontsize=25, fontweight='bold')
plt.legend(loc=4,frameon=False)
plt.grid()
plt.show()
# Please check out amorphous_silicon_Tersoff_LAMMPS in the example folder for more detail.