User Tools

Site Tools


absorbtion_spectra_by_fermi_golden_rule

This is an old revision of the document!


Absorbtion spectra by Fermi golden Rule

There is a possibility to evaluate rught estimation of UV-VIS absorption spectra by Fermi golden rule between Kohn-Sham orbitals. For given range range of photon energy are selected all pairs of occupied and unoccupied gorund-state Kohn-Sham orbital with appropriate energy difference. This two orbitals are than projected onto real space grid and transition dipole moment is computed. Oscilator strength (or probability of transition) is than proportional to squere of length of transition dipole moment. Absorbtion spectra can be than plotted by summing up contributions of oscillator strength for occupied-unoccupied orbital pairs with energy difference falling into particular energy bins (like weighted histogram).

Note that this method ( Fermi golden rule between ground state single particle states obtained from LDA-DFT ) is very crude approximation of true optical transitions. No beyond-DFT nor beyond-Hratree-Fock correction (such as RPA, TD-DFT, GW, BSE … ) are used. For this reason transition energy ( positioning of peak ) is mostly not realistic. The method is, however, useful to explore dependence of electronic transition probability on symmetry breaking and changes in shape and localization of electronic states.

Availability

this is implemented only in Prokop's personal versions of Fireball located in

/data/home/hapala/Fireball_dev/src_1.0-BSfinal
/data/home/hapala/Fireball_dev/progs_Jellium_mod

Usage

The computation of transitions is activated by keyword iwrtexcit = 1 in firebal.in. However, currently it is also necessary to activate iwrtewf = 1 in order to initialize grid, even thought .xsf files written by iwrtewf option are not used at all and has nothing incommion with the computatuion of transitions ( this should be corrected in future). So the typical fireball.in looks like:

&OPTION
basisfile  = answer.bas 
nstepf     = 1
icluster   = 1
ifixcharge = 1
&END
&OUTPUT
iwrtewf     = 1
iwrtexcit=1
&END
&MESH
iewform = 3
npbands = 1
pbands = 1
&END

Than excitations.optional file should be provided to

 0        ! write_out_transition_map ( default 0, moslty just for debuging )
 0.0      ! dEmin ... minimum energy difference ( mostly 0)
 3.5      ! dEmax ... maximim energy difference

if write_out_transition_map is set to '1', an '.xsf' file containing spatial distribution of transition dipole moment is outputed for each pair of orbitals. This may be usefull to identify parts of molecule which contribute most to the electronic transition. However, the number of transitions grow quadratically with energy range and each '.xsf' file needs typically more than 100MB of disk space. So in most of cases you would quickly consume HUGE amount of disk space ( like 1 terrabyte ) if you switch this option on.

It is also recommended to chose reasonable small energy interval (dEmin,dEmax) oterwise the computation can take very long time.

Outputs

Excitations.dat file is writen out from the calculation. The file contains several columns of following meaingin:

  1. j,j index of ocupied and unoccupied orbital
  2. elta_Eij energy diffeence
  3. |D| length of transition dipole moment
  4. D_x,D_y,D_z cartesian components of transition dipole moment
  5. Ei,Ej energy eigenvalues corresponging to the Khon-Sham orbitals
   i      j         delta_Eij               |D|                 D_x                 D_y                   D_z               Ei                     Ej
  749    879        3.4862376070        0.9824594517        0.2759973643       -0.9336816482        0.1314937604       -3.0970434024       -6.5832810094
  750    879        3.4791174668        0.9753364513       -0.4413174043        0.2404091902        0.8358968616       -3.0970434024       -6.5761608692
  751    879        3.4582785159        0.5954896260        0.3594301771       -0.1688715641       -0.4437344220       -3.0970434024       -6.5553219183
  752    879        3.4453548553        0.8057311791       -0.1487142116       -0.4625862731       -0.6427291469       -3.0970434024       -6.5423982578
  753    879        3.4410779875        0.7684644226        0.0080712081       -0.6875508088       -0.3431418216       -3.0970434024       -6.5381213900
  754    879        3.4257564875        0.4202213606        0.0892087920       -0.1778280360       -0.3701418282       -3.0970434024       -6.5227998899

Making spectra

For transformation of pairwise table of transition dipole moment into optical spectra for plotting can be used for example following bash/awk scripts.

First script spectranice.sh accumulate (sum up) pairwise contributions of transition dipole moment into bins of given size and then blur them with gaussian filter (to make nice smooth spectra). The most importaint parameter is step size. For example with step size 0.001 eV:

spectranice.sh 0 3.5 0.001

spectranice.sh code:

min=$1
max=$2
step=$3
awk -v min=$min -v max=$max -v step=$step '
BEGIN{for(i=0;min+i*step<max;i++){hist[i]=0}; minx=10000000; maxx=-10000000  }
(NR>0){
x=$3; y=$4; yy=y*y;
ei=int((x-min)/step)
if(x>maxx){maxx=x}
if(x<minx){minx=x}
hist[ei-5]+=yy*0.00098
hist[ei-4]+=yy*0.00977
hist[ei-3]+=yy*0.04395
hist[ei-2]+=yy*0.11719
hist[ei-1]+=yy*0.20508
hist[ei+0]+=yy*0.24609
hist[ei+1]+=yy*0.20508
hist[ei+2]+=yy*0.11719
hist[ei+3]+=yy*0.04395
hist[ei+4]+=yy*0.00977
hist[ei+5]+=yy*0.00098
}
END{
#print maxx,minx
for(i=0;min+i*step<max;i++){
if(((min+i*step)>minx) && ((min+i*step)<maxx) ) {print min+i*step,hist[i]}}
}
' Excitations.dat > spectranice.dat

in some cases you are not interested in the sum of all cotributions ( it would correspond to absorbtion spectra) but rather in the transition with highest oscillator strength in given interval of energy (energy bin). For this purpose script 'spectramax.sh can be used.

spectramax.sh 0.001
#! /bin/bash
step=$1
min=`awk '
BEGIN{Emin=100000}
(NR>0){if($3<Emin)Emin=$3}
END{print Emin}
' Excitations.dat`
max=`awk '
BEGIN{Emax=0}
(NR>0){if($3>Emax)Emax=$3}
END{print Emax}
' Excitations.dat`
echo $min $max
awk -v min=$min -v max=$max -v step=$step '
BEGIN{for(i=0;min+i*step<max;i++){hist[i]=0;num[i]=0;}}
(NR>0){
ei=int(($3-min)/step)
if(hist[ei]<($4*$4)){hist[ei]=($4*$4)}
} 
END{for(i=0;min+i*step<max;i++){print min+i*step,hist[i]}}
' Excitations.dat > spectramax.dat
absorbtion_spectra_by_fermi_golden_rule.1406115797.txt.gz · Last modified: 2014/07/23 13:43 (external edit)