Shift current: Difference between revisions

From The Yambo Project
Jump to navigation Jump to search
(Blanked the page)
Tag: Blanking
 
(25 intermediate revisions by the same user not shown)
Line 1: Line 1:
[[File:Shift current.png|right|200px | Shift_current]]
== Introduction ==
In this tutorial we will show how to calculate Shift Current in bulk materials.<br>
We suppose you are already familiar with the non-linear response using the Yambo code.
If it is not the case please study the previous tutorials: <br>[[Linear response using Dynamical Berry Phase]] and [[Real time approach to non-linear response (SHG)]].<br>
This tutorial was created with the help of [https://www.linkedin.com/in/yuncheng-mao-8aba6279/ Yuncheng Mao].
== Setup calculations ==


In this tutorial we will take as example the two dimensional hBN.<br>
DFT wave-functions and inputs can be downloaded here: [http://media.yambo-code.eu/educational/tutorials/files/hBN-2D-RT.tar.gz hBN-2D-RT.tar.gz].<br>
First of all run the setup, then remove symmetries along the '''y''' direction, as explained in the tutorial above.
== Real-time setup and calculations ==
In order to generate input file for shift current you do: <span style="color:blue">yambo_nl -u n -V par -F input.in</span> <br>
nloptics                        # [R] Non-linear spectroscopy
NLogCPUs=0                      # [PARALLEL] Live-timing CPU`s (0 for all)
PAR_def_mode= "balanced"        # [PARALLEL] Default distribution mode ("balanced"/"memory"/"workload"/"KQmemory")
NL_CPU= " <span style="color:red">8 1</span>"                      # [PARALLEL] CPUs for each role
NL_ROLEs= " <span style="color:red">w k</span>"                    # [PARALLEL] CPUs roles (w,k)
DIP_CPU= ""                      # [PARALLEL] CPUs for each role
DIP_ROLEs= ""                    # [PARALLEL] CPUs roles (k,c,v)
OSCLL_CPU= ""                    # [PARALLEL] CPUs for each role
OSCLL_ROLEs= ""                  # [PARALLEL] CPUs roles (k,b)
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                    # [OPENMP/NL] Number of threads for nl-optics
OSCLL_Threads=0                  # [OPENMP/X] Number of threads for Oscillators
% NLBands
  <span style="color:red">17 |  24 | </span>                          # [NL] Bands range
%
NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLtime=-1.000000          fs    # [NL] Simulation Time
NLintegrator= "INVINT"          # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
NLCorrelation= "IPA"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX/LSEX/LHF")
NLLrcAlpha= 0.000000            # [NL] Long Range Correction
% NLEnRange
<span style="color:red"> 1.000000 | 8.000000 </span> |        eV    # [NL] Energy range (for loop on frequencies NLEnSteps/=0
%
NLEnSteps=  <span style="color:red">24    </span>              # [NL] Energy steps for the loop on frequencies
% NLrotaxis
  0.000000 | 0.000000 | 0.000000 |        # [NL] Rotation axis (for the loop on angles NLAngSteps/=0)
%
NLAngSteps=0                    # [NL] Angular steps (if NLAngSteps/=0 field versor will be ignored)
NLDamping= 0.200000        eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative RADLifeTime=Phase_LifeTime)
<span style="color:red">EvalCurrent</span>                    # [NL] Evaluate the current
#FrPolPerdic                  # [DIP] Force periodicity of polarization respect to the external field
% Field1_Freq
  0.100000 | 0.100000 |        eV    # [RT Field1] Frequency
%
Field1_NFreqs= 1                # [RT Field1] Frequency
Field1_Int=  1000.00      kWLm2 # [RT Field1] Intensity
Field1_Width= 0.000000    fs    # [RT Field1] Width
Field1_kind= "<span style="color:red">SOFTSIN</span>"          # [RT Field1] Kind(SIN|SOFTSIN| see more on src/modules/mod_fields.F)
Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
% Field1_Dir
<span style="color:red"> 0.000000 | 1.000000 | 0.000000 |  </span>    # [RT Field1] Versor
%
Field1_Tstart= 0.010000    fs    # [RT Field1] Initial Time
Notice that in this input we turned one the evaluation of current  <span style="color:red">EvalCurrent</span> and force the parallelization on the frequencies, <br><code>NL_CPU= " <span style="color:red">8 1</span></code>, <code> NL_ROLEs= " <span style="color:red">w k</span> </code>, that is much more efficient than the one on k-points.<br><br>
Then you can run simulation by doing: <span style="color:blue">yambo_nl -F input.in</span>
== Analysis of the results ==
In order to analyze the result we use the branch '''shift_current''' of this Yambopy: [https://github.com/attacc/yambopy/tree/shift_current https://github.com/attacc/yambopy/tree/shift_current],<br> and the following script for the post-processing:
from yambopy import *
from yambopy.nl.harmonic_analysis import Harmonic_Analysis
from yambopy.units import fs2aut
X_order=4
NLDB=YamboNLDB()
Harmonic_Analysis(NLDB,X_order=X_order,prn_Peff=True)
Then you can plot the file <code>o.YamboPy-Sigma_probe_order_0</code> that contains the DC part of the current. In the figure below I plot the 5th column<br>
that corresponds to the real part of <math>\sigma^{yyy}</math>:
[[File:Sigma2.png|center| 600px | Shift current]]
Then we used this script to plot the result:
import matplotlib.pyplot as plt
import numpy as np
fig = plt.figure(figsize=(10,10))
ax = fig.add_subplot()
plt.title("$|\sigma^{(2)}(0,\omega_1-\omega_1$)|",fontsize=20)
fontsize2=14
XHI0 = np.genfromtxt('o.YamboPy-Sigma_probe_order_0',comments="#")
plt.tick_params(\
    axis='x',          # changes apply to the x-axis
    which='both',      # both major and minor ticks are affected
    bottom='off',      # ticks along the bottom edge are off
    top='off',        # ticks along the top edge are off
    labelbottom='off')
plt.xlim(1.0,7.5)
plt.xlabel('Energy [eV]',fontsize=20)
plt.tick_params(axis='y', which='major', labelsize=14)
# plt.plot(XHI0[:,0], XHI0[:,2], linewidth=2.5, linestyle="-",label='$\sigma^{(xyy)}$')
plt.plot(XHI0[:,0], XHI0[:,4], linewidth=2.5, linestyle="-",label='Re $\sigma^{(yyy)}$')
plt.axhline(y=0.0, color='black', linestyle='-')
plt.legend(loc="upper right",fontsize=fontsize2)
plt.xticks(fontsize=fontsize2)
plt.savefig('sigma2.png', format='png')
plt.show()

Latest revision as of 18:38, 27 October 2025