Pump and Probe: Difference between revisions

From The Yambo Project
Jump to navigation Jump to search
No edit summary
 
(87 intermediate revisions by the same user not shown)
Line 1: Line 1:
'''This tutorial works only with Yambo version > 5.5, that will be released soon.'''


In this tutorial we will show you how to setup two external fields in yambo_nl, to perform pump and probe simulation. <br>
In this tutorial we will show you how to setup two external fields in '''yambo_nl''', to perform pump and probe simulation. <ref name="pnp">[https://arxiv.org/abs/2211.12241 Exciton - Exciton transitions involving strongly bound excitons: an ab-initio approach], D. Sangalli, M. D'Alessandro and C. Attaccalite, Physical Review B '''107''' (20), 205203 (2023) </ref><br>
The input simulation parameters of these tutorial require large computational power,
This tutorial is at independent particle level, but the same kind of simulation can be repeated including local field effects,  
we advice you to run them in parallel, or to use simplified parameters for example less bands, k-points and so on.
quasi-particle correction or electron-hole interaction, see for example the tutorials:<br>
This tutorial supposes that you are already familiar with real-time simulation with Yambo,
<br>
[[Correlation effects in the non-linear response]] <br>
[[Real time Bethe-Salpeter Equation (TDSE)]]<br>
<br>
This tutorial supposes that you are already familiar with real-time simulation with Yambo<ref name="dberry">[https://arxiv.org/abs/1609.09639 Non-linear response in extended systems: a real-time approach], C. Attaccalite </ref>
if it is not the case please check these tutorials: [[Tutorials#Non_linear_response]]
if it is not the case please check these tutorials: [[Tutorials#Non_linear_response]]


We start from the h-BN monolayer with an in place lattice constant a=2.5 Å and a box large 30 a.u. in the z-direction.
We start from the h-BN monolayer with an in place lattice constant a=2.5 Å and a box large 33 a.u. in the z-direction.
Standard DFT input for hBN monolayer for ABINIT or QuantumEspresso can be found here [[Tutorials]]. <br>
Standard DFT input for hBN monolayer for ABINIT or QuantumEspresso can be found here [[Tutorials]]. <br>
We used  a 24x24x1 k-points grid and 100 bands in the non-self consistent calculation.
We used  a 9x9x1 k-points grid and 8 bands in the non-self consistent calculation and 40 Ry cutoff for the wave-function.
In order to have faster calculation we advice you to decrease the k-points grid to 12x12x1, 
results will be very close to the converged ones.


In our example we choose direction [0,1,0] for the pump and the probe. We generate the ypp.in to remove symmetries with the command <span style="color:blue">ypp -y</span> and we modify it as:
In our example we choose direction [0,1,0] for the pump and [1,0,0] for the probe.  
We generate the ypp.in to remove symmetries with the command <span style="color:blue">ypp -y</span> and we modify it as:


  fixsyms                          # [R] Remove symmetries not consistent with an external perturbation
  fixsyms                          # [R] Remove symmetries not consistent with an external perturbation
Line 20: Line 22:
  %
  %
  % Efield2
  % Efield2
  <span style="color:red"> 0.000000 | 1.000000 | 0.000000 | </span>      # Additional external Electric Field
  <span style="color:red"> 1.000000 | 0.000000 | 0.000000 | </span>      # Additional external Electric Field
  %
  %
  BField= 0.000000          T    # [MAG] Magnetic field modulus
  BField= 0.000000          T    # [MAG] Magnetic field modulus
Line 29: Line 31:
  #RmSpaceInv                    # Remove Spatial Inversion
  #RmSpaceInv                    # Remove Spatial Inversion


then we go in the FixSymm directory and run again the setup. Now we calculate collisions,
then we go in the FixSymm directory and run again the setup. <br>
since we want to include excitonic effects in the calculations. <br>Using the command: <span style="color:blue">yambo_nl -v h+sex+cvonly -e -X s -r</span> we get the input file:
We will perform three different real-time calculations: <br>
* the probe field alone;
* the pump field alone;
* the pump and probe configuration.


collisions                      # [R] Collisions
== Probe only ==
em1s                            # [R][Xs] Statically Screened Interaction
dipoles                          # [R] Oscillator strenghts (or dipoles)
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
X_Threads=0                      # [OPENMP/X] Number of threads for response functions
RT_Threads=0                    # [OPENMP/RT] Number of threads for real-time
RandQpts= <span style="color:red">30000000</span>              # [RIM] Number of random q-points in the BZ
RandGvec= 1                RL    # [RIM] Coulomb interaction RS components
CUTGeo= " <span style="color:red">box Z</span>"                  # [CUT] Coulomb Cutoff geometry: box/cylinder/sphere/ws/slab X/Y/Z/XY..
% CUTBox
  0.00000 |  0.00000 |  <span style="color:red">29.50000</span> |        # [CUT] [au] Box sides
%
Chimod= "HARTREE"                # [X] IP/Hartree/ALDA/LRC/PF/BSfxc
% BndsRnXs
    1 | 100 |                        # [Xs] Polarization function bands
%
NGsBlkXs=  <span style="color:red"> 10000 mHa  </span>  # [Xs] Response block size
% LongDrXs
  0.000000 | 1.000000 | 0.000000 |        # [Xs] [cc] Electric Field
%
XTermKind= <span style="color:red">"BG" </span>              # [X] X terminator ("none","BG" Bruneval-Gonze)
% COLLBands
    <span style="color:red">3 | 6 |  </span>                        # [COLL] Bands for the collisions
%
HXC_Potential= "SEX+HARTREE+CVONLY" # [SC] SC HXC Potential
HARRLvcs= <span style="color:red">10000 mHa  </span> # [HA] Hartree    RL components
EXXRLvcs= <span style="color:red">10000 mHa  </span> # [XX] Exchange    RL components
CORRLvcs= <span style="color:red">10000 mHa  </span> # [GW] Correlation RL components


Notice that in the previous input file we added the flag CVONLY that force Yambo to calculate collisions integrals only between valence and conductions bands. This approximation is usually safe and speed up calculation, but it could fail in particular non-linear phenomena. <br>
We can generate the input file for a generic Pump and Probe calculations with the command: <span style="color:blue">yambo_nl -u p -F yambo.in_probe</span>:
Now that we calculated the COLLISIONS we will perform three different real-time calculations: 1) the probe field alone; 2) the pump field alone; 3) the pump and probe configuration.<br>
We can generate the input file for a generic Pump and Probe calculations with the command: <span style="color:blue">yambo_nl -u p -V qp -F yambo.in_probe</span>:
   
   
  nloptics                        # [R] Non-linear spectroscopy
  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")
DIP_CPU= ""                      # [PARALLEL] CPUs for each role
DIP_ROLEs= ""                    # [PARALLEL] CPUs roles (k,c,v)
  DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
  DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
  NL_Threads=0                    # [OPENMP/NL] Number of threads for nl-optics
  NL_Threads=0                    # [OPENMP/NL] Number of threads for nl-optics
  % NLBands
  % NLBands
   3 |  6 |                           # [NL] Bands range
   <span style="color:red">  3 |  6 |   </span>                      # [NL] Bands range
  %
  %
  NLverbosity= "high"              # [NL] Verbosity level (low | high)
  NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLstep= <span style="color:red">5.00000  </span>        as    # [NL] Time step length
  NLtime= <span style="color:red"> 80.000000</span>           fs    # [NL] Simulation Time
  NLtime= <span style="color:red">80.000000 </span>         fs    # [NL] Simulation Time
  NLintegrator= " <span style="color:red"> CRANKNIC</span>"           # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
  NLintegrator= "<span style="color:red">CRANKNIC</span>"         # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
  NLCorrelation= "IPA"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
  NLCorrelation= "<span style="color:red">SEX</span>"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
  NLLrcAlpha= 0.000000            # [NL] Long Range Correction
  NLLrcAlpha= 0.000000            # [NL] Long Range Correction
  NLDamping= 0.000000       eV    # [NL] Damping (or dephasing)
  NLDamping= <span style="color:red">0.000000 </span>      eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
  #EvalCurrent                  # [NL] Evaluate the current
  #EvalCurrent                  # [NL] Evaluate the current
  HARRLvcs= 913            RL    # [HA] Hartree    RL components
#FrPolPerdic                  # [DIP] Force periodicity of polarization respect to the external field
  EXXRLvcs= 913            RL    # [XX] Exchange    RL components
  HARRLvcs= 21817            RL    # [HA] Hartree    RL components
% GfnQP_E
  EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
  <span style="color:red">3.594325</span> | 1.000000 | 1.000000 |        # [EXTQP G] E parameters  (c/v) eV|adim|adim
%
  % Field1_Freq
  % Field1_Freq
   0.100000 | 0.100000 |        eV    # [RT Field1] Frequency
   0.100000 | 0.100000 |        eV    # [RT Field1] Frequency
  %
  %
  Field1_Int=  <span style="color:red">1000.00</span>       kWLm2 # [RT Field1] Intensity
Field1_NFreqs= 1                # [RT Field1] Frequency
  Field1_Int=  1000.00      kWLm2 # [RT Field1] Intensity
  Field1_Width= 0.000000    fs    # [RT Field1] Width
  Field1_Width= 0.000000    fs    # [RT Field1] Width
  Field1_kind= "<span style="color:red">DELTA</span>"             # [RT Field1] Kind(SIN|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field1_kind= " <span style="color:red"> DELTA</span>"           # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
  Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
  % Field1_Dir
  % Field1_Dir
  <span style="color:red">0.000000 | 1.000000 | 0.000000 | </span>       # [RT Field1] Versor
  <span style="color:red"> 1.000000 | 0.000000 | 0.000000 | </span>     # [RT Field1] Versor
  %
  %
% Field1_Dir_circ
  Field1_Tstart= <span style="color:red"> 5.00000 </span> fs    # [RT Field1] Initial Time
  0.000000 | 1.000000 | 0.000000 |        # [RT Field1] Versor_circ
[......]
%
  Field1_Tstart= <span style="color:red">5.00000</span>   fs    # [RT Field1] Initial Time


Notice that we introduced a scissor operator in such a way to have the main excitonic peak at 6.1 eV in agreement with experimental measurements. In this input we set only the first field that is the probe, and we use a DELTA function as field in such a way to probe all the spectra. Then we can analyze the spectra with <span style="color:blue">ypp -u</span>:
In this input we set only the first field that is the probe, we use a DELTA function as field in such a way to probe all the frequencies, and set the starting point of the delta at t=5 fs. <br>Run this job in a separate folder with the command  <span style="color:blue">yambo_nl -F yambo.in_probe -J PROBE</span>. Then we can analyze the spectra with <span style="color:blue">ypp -u -J PROBE </span>:


  nonlinear                        # [R] Non-linear response analysis
  nonlinear                        # [R] Non-linear response analysis
Line 114: Line 83:
   0.00000 | 20.00000 |        eV    # Energy range
   0.00000 | 20.00000 |        eV    # Energy range
  %
  %
  DampMode= "LORENTZIAN"          # Damping type ( NONE | LORENTZIAN | GAUSSIAN )
  DampMode= "<span style="color:red">LORENTZIAN</span>"          # Damping type ( NONE | LORENTZIAN | GAUSSIAN )
  DampFactor= 0.100000      eV    # Damping parameter
  DampFactor= <span style="color:red">0.100000</span>       eV    # Damping parameter
  PumpPATH= "none"          # Path of the simulation with the Pump only
  PumpPATH= "none"          # Path of the simulation with the Pump only


and get the dielectric response. Hereafter we plot the polarization along the y-direction and the corresponding dielectric constant:
and get the dielectric response. Notice that we introduce a Lorentzian broadening the post-processing. <br>
Hereafter we plot the damped polarization along the y-direction and the corresponding dielectric constant:
 
[[File:Eps probe.png|800px|center | Probe only]]
 
the arrow indicates the position where we will put the Pump in the next simulation.


[[File:Pol eps small.png|800px|center|Yambo tutorial image]]
==Pump only==


Now we start from the previous input <span style="color:green">cp yambo.in_probe yambo.in_pump</span> and <br> modify it
Now we run a simulation with the pump field only. We generate the input in the same way <span style="color:blue">yambo_nl -u p -F yambo.in_pump</span>. <br>
in order perform a calculation with the Pump field only<span style="color:blue">yambo_nl -u p -V nl -F yambo.in_pump</span>. These result will be our reference field in the final analysis:
These result will be our reference field in the final analysis. Here the input file for the pump:
   
   
  nloptics                        # [R] Non-linear spectroscopy
  nloptics                        # [R] Non-linear spectroscopy
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                    # [OPENMP/NL] Number of threads for nl-optics
  % NLBands
  % NLBands
  3 |  6 |                          # [NL] Bands range  
    <span style="color:red">3 |  6 |</span>                           # [NL] Bands range
  %
  %
  NLverbosity= "high"              # [NL] Verbosity level (low | high)
  NLverbosity= "high"              # [NL] Verbosity level (low | high)
  NLstep= 5.00000          as    # [NL] Time step length
  NLtime<span style="color:red">150.0000</span>           fs    # [NL] Simulation Time
  NLtime= 150.000000           fs    # [NL] Simulation Time
  NLintegrator= "CRANKNIC"        # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
  NLintegrator= "CRANKNIC"        # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
  NLCorrelation= "SEX"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
  NLCorrelation= "IPA"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
  NLLrcAlpha= 0.000000            # [NL] Long Range Correction
  NLLrcAlpha= 0.000000            # [NL] Long Range Correction
  NLDamping= 0.000000        eV    # [NL] Damping (or dephasing)
  NLDamping= <span style="color:red">0.000000</span>       eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
  #EvalCurrent                  # [NL] Evaluate the current
  #EvalCurrent                  # [NL] Evaluate the current
  <span style="color:red">NoPperiodic</span>                    # [DIP] Do not enforce periodicity for the Polarization
  #FrPolPerdic                  # [DIP] Force periodicity of polarization respect to the external field
  HARRLvcs= 913            RL    # [HA] Hartree    RL components
  HARRLvcs= 21817            RL    # [HA] Hartree    RL components
  EXXRLvcs= 913            RL    # [XX] Exchange    RL components
  EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
  % GfnQP_E
  % Field1_Freq
   3.594325 | 1.000000 | 1.000000 |        # [EXTQP G] E parameters (c/v) eV|adim|adim
   0.100000 | 0.100000 |         eV    # [RT Field1] Frequency
%
Field1_NFreqs= 1                 # [RT Field1] Frequency
Field1_Int=  10000.00      kWLm2 # [RT Field1] Intensity
Field1_Width= 0.000000    fs    # [RT Field1] Width
Field1_kind= "SOFTSIN"          # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
% Field1_Dir
  0.000000 | 0.000000 | 0.000000 |        # [RT Field1] Versor
  %
Field1_Tstart= 0.010000    fs    # [RT Field1] Initial Time
% Field2_Freq
  <span style="color:red"> 5.42 | 5.42 </span>|        eV    # [RT Field2] Frequency
%
Field2_NFreqs= 1                # [RT Field2] Frequency
Field2_Int=  <span style="color:red">1000.00 </span>      kWLm2 # [RT Field2] Intensity
Field2_Width= <span style="color:red">10.00000</span>    fs    # [RT Field2] Width
Field2_kind= "<span style="color:red">QSSIN</span>"            # [RT Field2] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field2_pol= "linear"            # [RT Field2] Pol(linear|circular)
% Field2_Dir
<span style="color:red"> 0.000000 | 1.000000 | 0.000000 |  </span>    # [RT Field2] Versor
%
Field2_Tstart= 0.010000    fs    # [RT Field2] Initial Time
[...]
In this simulation we used as Pump field QSSIN that is sinus with a Guassian envelope[https://demonstrations.wolfram.com/SineGaussianSignals/ sine-Gaussian wave]. For a list of all possible field in Yambo see the file  <span style="color:green"> src/modules/mod_fields.F</span>. We centred the field frequency at 5.42 eV that is the energy of the highest peak in the previous figure, see arrow there.
Run the simulation with the command <span style="color:blue">yambo_nl -F yambo.in_pump -J PUMP </span>.
The pump field, column 6th of the file <code>o-PUMP.external_potential_F1</code>, is shown in the figure below:
 
[[File:Pump only.png|800px|center|Pump field]]
 
== Pump and probe==
 
Now we copy the pump in file in a new file <span style="color:green"> cp yambo.in_pump yambo.in_pump_and_probe</span> and turn on  the probe field:
 
nloptics                        # [R] Non-linear spectroscopy
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                    # [OPENMP/NL] Number of threads for nl-optics
% NLBands
  3 | 6 |                           # [NL] Bands range
  %
  %
NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLtime=  <span style="color:red">150.0000 </span>          fs    # [NL] Simulation Time
NLintegrator= "CRANKNIC"        # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
NLCorrelation= "IPA"            # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
NLLrcAlpha= 0.000000            # [NL] Long Range Correction
NLDamping=  <span style="color:red">0.000000 </span>      eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
#EvalCurrent                  # [NL] Evaluate the current
#FrPolPerdic                  # [DIP] Force periodicity of polarization respect to the external field
HARRLvcs= 21817            RL    # [HA] Hartree    RL components
EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
  % Field1_Freq
  % Field1_Freq
   0.100000 | 0.100000 |        eV    # [RT Field1] Frequency
   0.100000 | 0.100000 |        eV    # [RT Field1] Frequency
  %
  %
  Field1_Int=  1000.00       kWLm2 # [RT Field1] Intensity
Field1_NFreqs= 1                # [RT Field1] Frequency
  Field1_Int=  <span style="color:red">1.00 </span>    kWLm2 # [RT Field1] Intensity
  Field1_Width= 0.000000    fs    # [RT Field1] Width
  Field1_Width= 0.000000    fs    # [RT Field1] Width
  Field1_kind= "DELTA"             # [RT Field1] Kind(SIN|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field1_kind= "<span style="color:red">DELTA</span>"           # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
  Field1_pol= "linear"            # [RT Field1] Pol(linear|circular)
  % Field1_Dir
  % Field1_Dir
  0.000000 | 1.000000 | 0.000000 |        # [RT Field1] Versor
<span style="color:red"> 1.000000 | 0.000000 | 0.000000 | </span>      # [RT Field1] Versor
%
% Field1_Dir_circ
  0.000000 | 0.000000 | 0.000000 |       # [RT Field1] Versor_circ
  %
  %
  Field1_Tstart= <span style="color:red">165.00000    fs</span>    # [RT Field1] Initial Time
  Field1_Tstart= <span style="color:red">90.0000  </span> fs   # [RT Field1] Initial Time
  % Field2_Freq
  % Field2_Freq
  <span style="color:red">6.1</span> | 0.00000 |        eV    # [RT Field2] Frequency
  <span style="color:red">5.42 | 5.42 |  </span>       eV    # [RT Field2] Frequency
  %
  %
  Field2_Int= <span style="color:red">1.E4</span>     kWLm2 # [RT Field2] Intensity
Field2_NFreqs= 1                # [RT Field2] Frequency
  Field2_Width= <span style="color:red">10.000000    fs</span>    # [RT Field2] Width
  Field2_Int= <span style="color:red">1000.00</span>       kWLm2 # [RT Field2] Intensity
  Field2_kind= "<span style="color:red">QSSIN</span>"             # [RT Field2] Kind(SIN|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field2_Width= <span style="color:red">10.00000  </span>   fs   # [RT Field2] Width
  Field2_kind= "<span style="color:red">QSSIN</span>"             # [RT Field2] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
  Field2_pol= "linear"            # [RT Field2] Pol(linear|circular)
  Field2_pol= "linear"            # [RT Field2] Pol(linear|circular)
  % Field2_Dir
  % Field2_Dir
   0.000000 | 1.000000 | 0.000000 |        # [RT Field2] Versor
   <span style="color:red">0.000000 | 1.000000 | 0.000000 |</span>       # [RT Field2] Versor
%
Field2_Tstart= 0.010000    fs    # [RT Field2] Initial Time
 
'''Please note:'''  by convention in Yambo we consider the first field the probe and the second the pump. <br>
 
you run the simulation with the command <span style="color:blue">yambo_nl -F yambo.in_pump_and_probe -J PUMP_AND_PROBE </span> and then we analyze the result using ypp with the command: <span style="color:blue">ypp_nl -u -J PUMP_AND_PROBE </span>:
 
nonlinear                        # [R] Non-linear response analysis
Xorder= 1                        # Max order of the response/exc functions
% TimeRange
-1.000000 |-1.000000 |        fs    # Time-window where processing is done
  %
  %
  % Field2_Dir_circ
ETStpsRt= <span style="color:red">500 </span>                  # Total Energy steps
  0.000000 | 0.000000 | 0.000000 |       # [RT Field2] Versor_circ
  % EnRngeRt
  <span style="color:red">0.00000 | 3.9000 </span>|         eV    # Energy range
  %
  %
  Field2_Tstart= <span style="color:red">0.00000   fs</span>   # [RT Field2] Initial Time
  DampMode= "<span style="color:red">GAUSSIAN</span>"                # Damping type ( NONE | LORENTZIAN | GAUSSIAN )
DampFactor= <span style="color:red">0.085000  </span>    eV   # Damping parameter
PumpPATH= "<span style="color:red">./PUMP/</span>"                # Path of the simulation with the Pump only
 
Notice that in this case we are interested in the spectra below the 4 eV, that is was zero without the presence of the pump. We put in the input the pump path in such way to calculate polarization of the pump&probe minus the one of generate by the pump only <math>\Delta P(t) = P_{PnP}(t) - P_{pump}(t)</math>. Then we analyze the <math>\Delta P(t)</math> and get the new spectra:
 
[[File:P and p spectrum.png|800px | center|Pump and probe spectrum]]
 
The peak that you see in this spectrum, that were not present in when the pump was zero, correspond to the transition between the different excited states. The pump field populates the state at 5.42 eV and then the probe excites electrons from this state to higher states. The energy difference between the higher excited states and the peak at 5.42 define the position of the peaks in the pump and probe spectrum, while their intensity is dictated by the dipole matrix elements between the different excited states, see Ref. <ref name="pnp"></ref>.<br>
 
'''Additional comments on pump and probe calculations '''
* We used a Gaussian damping factor, this choice is motivated by various tests we have done, and it seems to us that in pump and probe this broadening gives better results.<br>
* If you plot all spectrum you will notice that the "equilibrium part" is much stronger than the one induced by the pump, this is normal because only a small fraction of electron are promoted in the excited state by the pump.
* All previous simulation can be run also in presence of excitons, see  Ref. <ref name="pnp"></ref> and group theory can be used to predict dark/bright transition in the pump and probe spectroscopy.<ref>[https://link.springer.com/book/10.1007/978-3-540-32899-5  Group Theory, Application to the Physics of Condensed Matter],  Mildred S. S. Dresselhaus, Gene Dresselhaus, Ado Jorio, Springer Editor</ref>
* In principle is it possible to include damping and radiative life-time in the simulation to study retardation effect between the pump and probe. This  approach  exact in the independent particle case, but introduces an error when correlation effects are present, see Refs.<ref>[https://arxiv.org/abs/2106.03592 Excitons and carriers in transient absorption and time-resolved ARPES spectroscopy: an abinitio approach], D. Sangalli, Phys. Rev. Materials '''5''', 083803 (2021)</ref> and <ref>[https://arxiv.org/abs/1705.04245 Theory of exciton-phonon coupling], Gabriel Antonius and Steven G. Louie Phys. Rev. B '''105''', 085111 (2022)</ref>for a discussion.
* Non-linear response (SHG, THG, ...) induced by the pump is not implemented yet... if you have a good idea how to do it, please let us know
 
== Analysis of the results using YamboPy ==
'''This part works only with Yambo 5.3 or the last version available on Github:''' [https://github.com/yambo-code/yambo https://github.com/yambo-code/yambo]
<br>
Results can be analyzed using YamboPy with the following script:
 
from yambopy import *
from yambopy.plot  import *
NLDB_P_and_P=YamboNLDB(calc='P_and_P')
NLDB_Pump  =YamboNLDB(calc='PUMP')
 
pol_pump  =NLDB_Pump.Polarization[0]
pol_p_and_p=NLDB_P_and_P.Polarization[0]
time=NLDB_P_and_P.IO_TIME_points
t_initial=NLDB_P_and_P.Efield[0]["initial_time"]
pol_damped=np.empty_like(pol_pump)
for i_d in range(3):
    pol_damped[i_d,:]=damp_it(pol_p_and_p[i_d,:]-pol_pump[i_d,:],time,t_initial,damp_type='GAUSSIAN',damp_factor=0.085/ha2ev)
Plot_Pol(time=time, pol=pol_damped, xlim=[0,55], save_file='polarization.pdf')
Linear_Response(time=time,pol=pol_damped,efield=NLDB_P_and_P.Efield[0],plot=False,plot_file='delta_eps.pdf',e_range=[0.0,3.7],n_freqs=500)


We set the starting time of the Probe (Field1) to 165 fs, a time larger then the simulation run in such a way to have only the Pump(Field2) field. As Pump field we use QSSIN that is an [https://demonstrations.wolfram.com/SineGaussianSignals/ sine-Gaussian wave], with a frequency 6.1 eV centered at the energy of the first exciton, and a width of 10 fs. Notice
== References ==

Latest revision as of 08:56, 2 September 2024

In this tutorial we will show you how to setup two external fields in yambo_nl, to perform pump and probe simulation. [1]
This tutorial is at independent particle level, but the same kind of simulation can be repeated including local field effects, quasi-particle correction or electron-hole interaction, see for example the tutorials:

Correlation effects in the non-linear response
Real time Bethe-Salpeter Equation (TDSE)

This tutorial supposes that you are already familiar with real-time simulation with Yambo[2] if it is not the case please check these tutorials: Tutorials#Non_linear_response

We start from the h-BN monolayer with an in place lattice constant a=2.5 Å and a box large 33 a.u. in the z-direction. Standard DFT input for hBN monolayer for ABINIT or QuantumEspresso can be found here Tutorials.
We used a 9x9x1 k-points grid and 8 bands in the non-self consistent calculation and 40 Ry cutoff for the wave-function.

In our example we choose direction [0,1,0] for the pump and [1,0,0] for the probe. We generate the ypp.in to remove symmetries with the command ypp -y and we modify it as:

fixsyms                          # [R] Remove symmetries not consistent with an external perturbation
% Efield1
 0.000000 | 1.000000 | 0.000000 |        # First external Electric Field
%
% Efield2
 1.000000 | 0.000000 | 0.000000 |        # Additional external Electric Field
%
BField= 0.000000           T     # [MAG] Magnetic field modulus
Bpsi= 0.000000             deg   # [MAG] Magnetic field psi angle [degree]
Btheta= 0.000000           deg   # [MAG] Magnetic field theta angle [degree]
#RmAllSymm                     # Remove all symmetries 
RmTimeRev                     # Remove Time Reversal
#RmSpaceInv                    # Remove Spatial Inversion

then we go in the FixSymm directory and run again the setup.
We will perform three different real-time calculations:

  • the probe field alone;
  • the pump field alone;
  • the pump and probe configuration.

Probe only

We can generate the input file for a generic Pump and Probe calculations with the command: yambo_nl -u p -F yambo.in_probe:

nloptics                         # [R] Non-linear spectroscopy
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                     # [OPENMP/NL] Number of threads for nl-optics
% NLBands
    3 |  6 |                           # [NL] Bands range
%
NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLtime=  80.000000           fs    # [NL] Simulation Time
NLintegrator= "  CRANKNIC"           # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
NLCorrelation= "IPA"             # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
NLLrcAlpha= 0.000000             # [NL] Long Range Correction
NLDamping= 0.000000        eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
#EvalCurrent                   # [NL] Evaluate the current
#FrPolPerdic                   # [DIP] Force periodicity of polarization respect to the external field
HARRLvcs= 21817            RL    # [HA] Hartree     RL components
EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
% 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= "  DELTA"           # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field1_pol= "linear"             # [RT Field1] Pol(linear|circular)
% Field1_Dir
   1.000000 | 0.000000 | 0.000000 |        # [RT Field1] Versor
%
Field1_Tstart=   5.00000    fs    # [RT Field1] Initial Time
[......]

In this input we set only the first field that is the probe, we use a DELTA function as field in such a way to probe all the frequencies, and set the starting point of the delta at t=5 fs.
Run this job in a separate folder with the command yambo_nl -F yambo.in_probe -J PROBE. Then we can analyze the spectra with ypp -u -J PROBE :

nonlinear                        # [R] Non-linear response analysis
Xorder= 1                        # Max order of the response/exc functions
% TimeRange
-1.000000 |-1.000000 |         fs    # Time-window where processing is done
%
ETStpsRt= 200                    # Total Energy steps
% EnRngeRt
  0.00000 | 20.00000 |         eV    # Energy range
%
DampMode= "LORENTZIAN"           # Damping type ( NONE | LORENTZIAN | GAUSSIAN )
DampFactor= 0.100000       eV    # Damping parameter
PumpPATH= "none"           # Path of the simulation with the Pump only

and get the dielectric response. Notice that we introduce a Lorentzian broadening the post-processing.
Hereafter we plot the damped polarization along the y-direction and the corresponding dielectric constant:

Probe only

the arrow indicates the position where we will put the Pump in the next simulation.

Pump only

Now we run a simulation with the pump field only. We generate the input in the same way yambo_nl -u p -F yambo.in_pump.
These result will be our reference field in the final analysis. Here the input file for the pump:

nloptics                         # [R] Non-linear spectroscopy
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                     # [OPENMP/NL] Number of threads for nl-optics
% NLBands
   3 |  6 |                           # [NL] Bands range
%
NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLtime=  150.0000           fs    # [NL] Simulation Time
NLintegrator= "CRANKNIC"         # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
NLCorrelation= "IPA"             # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
NLLrcAlpha= 0.000000             # [NL] Long Range Correction
NLDamping=  0.000000        eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
#EvalCurrent                   # [NL] Evaluate the current
#FrPolPerdic                   # [DIP] Force periodicity of polarization respect to the external field
HARRLvcs= 21817            RL    # [HA] Hartree     RL components
EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
% Field1_Freq
 0.100000 | 0.100000 |         eV    # [RT Field1] Frequency
%
Field1_NFreqs= 1                 # [RT Field1] Frequency
Field1_Int=  10000.00       kWLm2 # [RT Field1] Intensity
Field1_Width= 0.000000     fs    # [RT Field1] Width
Field1_kind= "SOFTSIN"           # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field1_pol= "linear"             # [RT Field1] Pol(linear|circular)
% Field1_Dir
 0.000000 | 0.000000 | 0.000000 |        # [RT Field1] Versor
%
Field1_Tstart= 0.010000    fs    # [RT Field1] Initial Time 
% Field2_Freq
  5.42 | 5.42 |         eV    # [RT Field2] Frequency
%
Field2_NFreqs= 1                 # [RT Field2] Frequency
Field2_Int=  1000.00       kWLm2 # [RT Field2] Intensity
Field2_Width= 10.00000     fs    # [RT Field2] Width
Field2_kind= "QSSIN"             # [RT Field2] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field2_pol= "linear"             # [RT Field2] Pol(linear|circular)
% Field2_Dir
 0.000000 | 1.000000 | 0.000000 |        # [RT Field2] Versor
%
Field2_Tstart= 0.010000    fs    # [RT Field2] Initial Time
[...]

In this simulation we used as Pump field QSSIN that is sinus with a Guassian envelopesine-Gaussian wave. For a list of all possible field in Yambo see the file src/modules/mod_fields.F. We centred the field frequency at 5.42 eV that is the energy of the highest peak in the previous figure, see arrow there. Run the simulation with the command yambo_nl -F yambo.in_pump -J PUMP . The pump field, column 6th of the file o-PUMP.external_potential_F1, is shown in the figure below:

Pump field

Pump and probe

Now we copy the pump in file in a new file cp yambo.in_pump yambo.in_pump_and_probe and turn on the probe field:

nloptics                         # [R] Non-linear spectroscopy
DIP_Threads=0                    # [OPENMP/X] Number of threads for dipoles
NL_Threads=0                     # [OPENMP/NL] Number of threads for nl-optics
% NLBands
  3 |  6 |                           # [NL] Bands range
%
NLverbosity= "high"              # [NL] Verbosity level (low | high)
NLtime=  150.0000           fs    # [NL] Simulation Time
NLintegrator= "CRANKNIC"         # [NL] Integrator ("EULEREXP/RK2/RK4/RK2EXP/HEUN/INVINT/CRANKNIC")
NLCorrelation= "IPA"             # [NL] Correlation ("IPA/HARTREE/TDDFT/LRC/LRW/JGM/SEX")
NLLrcAlpha= 0.000000             # [NL] Long Range Correction
NLDamping=  0.000000        eV    # [NL] Damping (or dephasing)
RADLifeTime=-1.000000      fs    # [RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)
#EvalCurrent                   # [NL] Evaluate the current
#FrPolPerdic                   # [DIP] Force periodicity of polarization respect to the external field
HARRLvcs= 21817            RL    # [HA] Hartree     RL components
EXXRLvcs= 21817            RL    # [XX] Exchange    RL components
% Field1_Freq
 0.100000 | 0.100000 |         eV    # [RT Field1] Frequency
%
Field1_NFreqs= 1                 # [RT Field1] Frequency
Field1_Int=  1.00       kWLm2 # [RT Field1] Intensity
Field1_Width= 0.000000     fs    # [RT Field1] Width
Field1_kind= "DELTA"           # [RT Field1] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field1_pol= "linear"             # [RT Field1] Pol(linear|circular)
% Field1_Dir
 1.000000 | 0.000000 | 0.000000 |        # [RT Field1] Versor
%
Field1_Tstart= 90.0000    fs    # [RT Field1] Initial Time
% Field2_Freq
  5.42 | 5.42 |         eV    # [RT Field2] Frequency
%
Field2_NFreqs= 1                 # [RT Field2] Frequency
Field2_Int=  1000.00       kWLm2 # [RT Field2] Intensity
Field2_Width=  10.00000     fs    # [RT Field2] Width
Field2_kind= "QSSIN"             # [RT Field2] Kind(SIN|COS|RES|ANTIRES|GAUSS|DELTA|QSSIN)
Field2_pol= "linear"             # [RT Field2] Pol(linear|circular)
% Field2_Dir
 0.000000 | 1.000000 | 0.000000 |        # [RT Field2] Versor
%
Field2_Tstart= 0.010000    fs    # [RT Field2] Initial Time

Please note: by convention in Yambo we consider the first field the probe and the second the pump.

you run the simulation with the command yambo_nl -F yambo.in_pump_and_probe -J PUMP_AND_PROBE and then we analyze the result using ypp with the command: ypp_nl -u -J PUMP_AND_PROBE :

nonlinear                        # [R] Non-linear response analysis 
Xorder= 1                        # Max order of the response/exc functions
% TimeRange
-1.000000 |-1.000000 |         fs    # Time-window where processing is done
%
ETStpsRt= 500                    # Total Energy steps
% EnRngeRt
  0.00000 | 3.9000 |         eV    # Energy range
%
DampMode= "GAUSSIAN"                 # Damping type ( NONE | LORENTZIAN | GAUSSIAN )
DampFactor= 0.085000       eV    # Damping parameter
PumpPATH= "./PUMP/"                 # Path of the simulation with the Pump only

Notice that in this case we are interested in the spectra below the 4 eV, that is was zero without the presence of the pump. We put in the input the pump path in such way to calculate polarization of the pump&probe minus the one of generate by the pump only [math]\displaystyle{ \Delta P(t) = P_{PnP}(t) - P_{pump}(t) }[/math]. Then we analyze the [math]\displaystyle{ \Delta P(t) }[/math] and get the new spectra:

Pump and probe spectrum

The peak that you see in this spectrum, that were not present in when the pump was zero, correspond to the transition between the different excited states. The pump field populates the state at 5.42 eV and then the probe excites electrons from this state to higher states. The energy difference between the higher excited states and the peak at 5.42 define the position of the peaks in the pump and probe spectrum, while their intensity is dictated by the dipole matrix elements between the different excited states, see Ref. [1].

Additional comments on pump and probe calculations

  • We used a Gaussian damping factor, this choice is motivated by various tests we have done, and it seems to us that in pump and probe this broadening gives better results.
  • If you plot all spectrum you will notice that the "equilibrium part" is much stronger than the one induced by the pump, this is normal because only a small fraction of electron are promoted in the excited state by the pump.
  • All previous simulation can be run also in presence of excitons, see Ref. [1] and group theory can be used to predict dark/bright transition in the pump and probe spectroscopy.[3]
  • In principle is it possible to include damping and radiative life-time in the simulation to study retardation effect between the pump and probe. This approach exact in the independent particle case, but introduces an error when correlation effects are present, see Refs.[4] and [5]for a discussion.
  • Non-linear response (SHG, THG, ...) induced by the pump is not implemented yet... if you have a good idea how to do it, please let us know

Analysis of the results using YamboPy

This part works only with Yambo 5.3 or the last version available on Github: https://github.com/yambo-code/yambo
Results can be analyzed using YamboPy with the following script:

from yambopy import *
from yambopy.plot  import *

NLDB_P_and_P=YamboNLDB(calc='P_and_P')
NLDB_Pump   =YamboNLDB(calc='PUMP')
 
pol_pump   =NLDB_Pump.Polarization[0]
pol_p_and_p=NLDB_P_and_P.Polarization[0]

time=NLDB_P_and_P.IO_TIME_points
t_initial=NLDB_P_and_P.Efield[0]["initial_time"]

pol_damped=np.empty_like(pol_pump)

for i_d in range(3):
    pol_damped[i_d,:]=damp_it(pol_p_and_p[i_d,:]-pol_pump[i_d,:],time,t_initial,damp_type='GAUSSIAN',damp_factor=0.085/ha2ev)

Plot_Pol(time=time, pol=pol_damped, xlim=[0,55], save_file='polarization.pdf')

Linear_Response(time=time,pol=pol_damped,efield=NLDB_P_and_P.Efield[0],plot=False,plot_file='delta_eps.pdf',e_range=[0.0,3.7],n_freqs=500)

References

  1. 1.0 1.1 1.2 Exciton - Exciton transitions involving strongly bound excitons: an ab-initio approach, D. Sangalli, M. D'Alessandro and C. Attaccalite, Physical Review B 107 (20), 205203 (2023)
  2. Non-linear response in extended systems: a real-time approach, C. Attaccalite
  3. Group Theory, Application to the Physics of Condensed Matter, Mildred S. S. Dresselhaus, Gene Dresselhaus, Ado Jorio, Springer Editor
  4. Excitons and carriers in transient absorption and time-resolved ARPES spectroscopy: an abinitio approach, D. Sangalli, Phys. Rev. Materials 5, 083803 (2021)
  5. Theory of exciton-phonon coupling, Gabriel Antonius and Steven G. Louie Phys. Rev. B 105, 085111 (2022)