This notebook can be run on Google Colab.

Open in Colab

In Colab, you can enable the GPU acceleration from Edit > Notebook Settings > Accelerator > GPU.

Amorphous silicon tutorial

Complie LAMMPS as shared-library of python (\sim 8 min)

[ ]:
# Sudo apt-get softwares
! apt-get update
! apt install -y cmake build-essential git ccache openmpi-bin libopenmpi-dev python3.10-venv
! echo "Sudo apt-get finishes!"

# Upgrade python packages
! pip install --upgrade pip
! pip install numpy torch scipy virtualenv psutil pandas tabulate mpi4py Cython
! echo "Python packages installation finishes!"

# Build lammps with cmake
%cd /content
!rm -rf lammps
! git clone https://github.com/lammps/lammps.git lammps
%cd /content/lammps
! rm -rf build
! mkdir build
%cd build
! cmake ../cmake -DLAMMPS_EXCEPTIONS=yes \
               -DBUILD_SHARED_LIBS=yes \
               -DMLIAP_ENABLE_PYTHON=yes \
               -DPKG_PYTHON=yes \
               -DPKG_MANYBODY=yes \
               -DPKG_KSPACE=yes \
               -DPKG_PHONON=yes \
               -DPYTHON_EXECUTABLE:FILEPATH=`which python`

# Complie lammps as share-libary of python
! make -j 2
! make install-python
! echo "LAMMPS compilation done!"

# Redirect back to main folder
%cd /content/

Remote fetch and install source code from Github

[ ]:
! git lfs install --skip-repo --skip-smudge
! pip install git+https://github.com/nanotheorygroup/kaldo

Fetch supplyment data remotely

[ ]:
# Remote fetch kaldo resources from drop box
! wget https://www.dropbox.com/s/bvw0qcxy397g25q/kaldo_resources.zip?dl=0
! mv kaldo_resources.zip?dl=0 kaldo_resources.zip
! unzip kaldo_resources.zip

# Unzip files
! unzip forcefields.zip
! unzip structure_a_si_512.zip

# Clean workspace
! rm -r forcefields.zip
! rm -r structure_a_si_512.zip
! rm -r kaldo_resources.zip

Thermal transport simulation for Amorphous-Silicon (a-Si)

[ ]:
from kaldo.forceconstants import ForceConstants

# Read amporphous silicon with eskm format
forceconstants = ForceConstants.from_folder(folder='structure_a_si_512', format='eskm')

Create phonons object

[5]:
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.05/4.135,
                   broadening_shape='triangle',
                   storage='numpy')

Access and visualize properties calculated during simulations

[ ]:
import matplotlib.pyplot as plt
plt.style.use('seaborn-poster')

# 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.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.show()

Calculate and visualize \kappa_{per \ mode} and \kappa_{cum}

[ ]:
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.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.