{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Temperature-dependent properties using the ZG method\n", "\n", "Here we calculate temperature-dependent properties such as band strictures using the ZG method.\n", "\n", "Below we define the main input file which defines all constants and the system we will investigate, in this case silicon.\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " \n", " \n", " -*#*- ...............- \n", " .+*= .+%*-=%%: .=#*- -===============-:.\n", " :*%=*%%- *%* #%* :+%+-%%+ .:. -=. :==-. \n", " -%S -%%*: :#%. -%%-. -##: #%* -=. :==- \n", " .. .%S: +%%%%*. :*%%%#= %%= -=. :==- \n", " :=#%%*- .#S- .. .%%= :*#*: -=. :==- \n", " -%S:.=#%%*==%# *%%=::=##-+%%. . .=-. :==- .= \n", " :%%- .-+++: -+##*=. =%S :-::==: .==- --.\n", " *%# #%+ -=--:. .----:. \n", " :%%- -%S. \n", " .-=*####SS%#########: -###########*+: +#######= =%SS####+::. \n", " =+#**%SSSSSSSSSSSSSSSS= =SSSSSSSSSSSSSSS= #SSSSSSS*. +SSSSSSSS%%%%- \n", " *%% =SSS.. .SSS= SSS=. .:+SSS+ -SSS:. .-::. .SSS+. #%* \n", " #%#. =SSS. *S# *S%: SSS= %SS%. .SSS= #SSS%. :SSS: =%#. \n", " *%%: =SSS#*#SSS- SSS= .+SSS+. %SS*.=SSSSS+ +SS%.:+%+ \n", " +%%:=SSSSSSSSS- SSSSSSSSSSSSS=. +SSS:SSS%SSS:%SS*-#%= \n", " ....#S==SSS:..SSS: =+- SSS%######+=. -SSS%SS%.#SS%SSS==%%=. \n", " .:+##%%#*- =SSS. ::. :SSS. SSS=. SSSSSS:..SSSSSS: :*%%%*=- \n", " #%+. -+#SSS*+++++++*SSS: .=+SSS#++++: %SSSS+. =SSSS%. :-+#%%+ \n", " #%* .SSSSSSSSSSSSSSSSSS: *SSSSSSSSSSS. +SSS%. %SSS*. . .%%= \n", " =%S :::::::::::::::::. .:::::::::. :::. :::. =+=-.#%+ \n", " -%S: =+===#S: \n", " ==*------------------------------=========+++++++++++++++++++++++========++-+## \n", " =+++++++++++*******++++++++++++++++========------------------==========+++++++- \n", "-- -- -- -- -- -- -- -- -- -- -- -- Structure Info -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "1\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si.upf\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si_r.upf\n", "pseudo found at pseudodojo : ONCVPSP-PBE-FR-PDv0.4/Si_r.upf\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "pseudopotential: Si_r.upf\n", "pseudopotential directory: '/mnt/storage/sabya/For_video/epwpy/notebooks_basic/pseudo/'\n", "prefix: si\n" ] } ], "source": [ "import numpy as np\n", "import math\n", "import os\n", "import sys\n", "import subprocess\n", "import EPWpy\n", "from EPWpy import EPWpy\n", "from EPWpy.plotting import plot_bands\n", "\n", "####Constants########\n", "nr=3.3\n", "hbar=6.6*10**-16\n", "c=3*10**10\n", "font=16\n", "# Supercell dimensions\n", "dim1=3\n", "dim2=3\n", "dim3=3\n", "\n", "#Define folder locations, prefix and number of cores\n", "folder='./'\n", "cores='16'\n", "prefix='si'\n", "flfrc='\\''+str(prefix)+'.fc\\''\n", "pseudo = os.getcwd() + '/pseudos'\n", "######Define the directory of installation##############\n", "QE = '/home1/07369/mzach/codes/q-e_dev_2024/bin'\n", "########################################################\n", "\n", "\n", "########Define the system used for the calculations###############\n", "silicon=EPWpy.EPWpy({'prefix':prefix,'restart_mode':'\\'from_scratch\\'','ibrav':2,'nat':2,'calculation':'\\'scf\\'',\n", " 'atomic_species':['Si'],'mass':[28.0855],\n", " 'atoms':['Si','Si'],'ntyp':1,'pseudo':['Si.upf'],'ecutwfc':'40','ecutrho':'160',\n", " 'celldm(1)':'10.262','verbosity':'high','type':'alat','pseudo_auto':True \n", " },env='mpirun')\n", "silicon.run_serial = True\n", "\n", "#######Printing any attribute######\n", "pseudopot=silicon.__dict__['pw_atomic_species']['pseudo'][0]\n", "print('pseudopotential:', silicon.__dict__['pw_atomic_species']['pseudo'][0])\n", "print('pseudopotential directory:', silicon.__dict__['pw_control']['pseudo_dir'])\n", "print('prefix:',prefix)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Self\\,Consistent\\,Field\\,(SCF) \\,Calculations}}$\n", "\n", "We first solve the Kohn-Sham equations to obtain the Kohn-Sham orbitals $\\phi_v(r)$, $r$ is the electronic position (generally a mesh grid), $R$ is the position of ions.\n", "\n", "$E[\\phi_v,R]=-\\frac{\\hbar^2}{2m}\\sum_v{\\int{\\phi_v^\\star(r)\\nabla^2\\phi_v(r)dr}+\\int{V_R(r)n(r)dr}+\\frac{e^2}{2}\\int{\\frac{n(r)n(r')}{|r-r'|}drdr'}+E_{xc}[n(r')]+\\sum_{I\\neq J}{\\frac{e^2}{2}\\frac{Z_IZ_J}{|R_I-R_J|}}}$\n", "\n", "We minimize $E(R)=min(E[\\phi_v,R])$\n", "\n", "Where, $\\Big(-\\frac{\\hbar^2}{2m}\\nabla^2+V_{KS}(r)\\Big)\\phi_v(r)=\\epsilon_v\\phi_v(r)$\n" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: scf -- -- -- -- -- -- -- -- -- -- -- \n", "Running scf |████████████████████████████████████████| in 1.4s (2.74/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Create SCF input file ##\n", "silicon.scf(electrons={'conv_thr':'1E-6',\n", " 'mixing_beta':'0.7'},\n", " kpoints={'kpoints':[[3,3,3]],\n", " 'kpoints_type':'automatic'},\n", " name='scf')\n", "## Prepare folders and copy necessary files ##\n", "silicon.prepare(1,type_run='scf')\n", "## Run calculation ##\n", "silicon.run(2)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "Calculation finished normally in si/scf/scf.out\n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "silicon.file = 'si/scf/scf.out'\n", "silicon.get_QE_status()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Phonon\\,Calculations}}$\n", "\n", "In order to obtain the temperature-dependent properties we need the electron-phonon interactions, which is the first derivative of the lattice potential.\n", "\n", "We also need phonon energies and displacement vectors. We will obtain all these quantities using the density-functional perturbation theory (DFPT)\n", "\n", "https://docs.epw-code.org/_downloads/b3f5899664a87fcdd6dcacc262e6f103/Mon.1.Giannozzi.pdf\n", "\n", "Here as an example we evaluate phonons in a $3\\times3\\times3$ $q$-grid." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: ph -- -- -- -- -- -- -- -- -- -- -- \n", "Running ph |████████████████████████████████████████| in 15.3s (0.09/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Create PH input file ##\n", "silicon.ph(phonons={'nq1':3,\n", " 'nq2':3,\n", " 'nq3':3,\n", " 'fildvscf':'\\' \\''})\n", "## Prepare folders and copy necessary files ##\n", "silicon.prepare(type_run='ph')\n", "## Run calculation ##\n", "silicon.run(4,type_run='ph')" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "Calculation finished normally in si/ph/ph.out\n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "silicon.file = 'si/ph/ph.out'\n", "silicon.get_QE_status()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Interatomic\\,Force\\,Constants\\,(IFCs)}}$\n", "\n", "In order to produce the ZG configurations we first need to evaluate the interatomic forces constants, as they will provide information on how to displace the atoms. This constants are obtain via second derivatives of the total potential energy, that is: $C_{\\kappa\\alpha\\rho\\kappa'\\alpha'\\rho'}=\\partial^{2}U/\\partial\\tau_{\\kappa\\alpha\\rho}\\tau_{\\kappa'\\alpha'\\rho'}$" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: q2r -- -- -- -- -- -- -- -- -- -- -- \n", "Running q2r |████████████████████████████████████████| in 1.2s (4.12/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Create q2r input file ##\n", "silicon.q2r(q2r={'fildyn':'\\'si.dyn\\'',\n", " 'flfrc':flfrc,\n", " 'zasr':'\\'crystal\\''})\n", "## Prepare folders and copy necessary files ##\n", "silicon.prepare(1,type_run='q2r')\n", "## Run calculation ##\n", "silicon.run(1,type_run='q2r')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm The\\,ZG\\,Displacements}}$\n", "\n", "In order to generate the set of atomic displacements (ZG displacements) that best incorporate the effect of electron-phonon coupling in ab initio nonperturbative calculations we employ the equation described below, as implemented in Quantum ESPRESSO by M. Zacharias and F. Giustino (Phys. Rev. Research 2, 013357 (2020)):\n", "$\\Delta\\tau^{ZG}_{pk}=\\sqrt{\\frac{M_{p}}{N_{p}M_{k}}}2\\sum_{q\\in B,v}S_{qv}Re[e^{iq.R_{p}}e_{k,v}(q)]\\sigma_{qv,T}$\n", "\n", "In this example we will create ZG configurations in a $3\\times3\\times3$ at 0 Kelvin. Our goal in this example is to determine the phonon-induced band gap renormalization." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "obtaining nscf and ph attributes\n", "-- -- -- -- -- -- -- -- -- -- -- Calculation: zg -- -- -- -- -- -- -- -- -- -- -- \n", "Running zg |████████████████████████████████████████| in 1.8s (1.44/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Create ZG input file ##\n", "silicon.reset()\n", "size=3 # Supercell size given by size X size X size\n", "T=0.0 # Temperature\n", "silicon.zg(zg={'T':T,\n", " 'flfrc':flfrc,\n", " 'dim1':dim1,\n", " 'dim2':dim2,\n", " 'dim3':dim3,\n", " 'error_thresh':'0.3',\n", " 'niters':'4000'})\n", "## Prepare folders and copy necessary files ##\n", "silicon.prepare(type_run='zg')\n", "## Run calculation ##\n", "silicon.run(2,type_run='zg')" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "Calculation finished normally in si/zg/zg.out\n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "silicon.file = 'si/zg/zg.out'\n", "silicon.get_QE_status()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Phonon-Assisted\\,Optical\\,Spectra}}$\n", "\n", "Next, we will calculate temperature-dependent optical spectra. We will evaluate the optical matrix elements in the independent-particle approximation. Our goal is to obtain the real and imaginary part of the dielectric function and then calculate the absorption coefficient for both the equilibrium and ZG structure.\n", "\n", "First, we will run a SCF and NSCF calculation for both systems:\n", "\n", "1) The NSCF should be done with the flag \"nosym=.true.\", as the code we will use to calculate absorption does not recognize symmetries (epsilon_Gaus.x).\n", "\n", "2) We need to include extra bands with the \"nbnd\" flag to account for transitions to conduction states that are higher in energy.\n", "\n", "3) The epsilon_Gaus.x code doesn't support the reduction of the k-points grid into the irreducible Brillouin zone. Therefore the PW run must be performed with a uniform k-points grid and all k-points weights must be equal to each other. In other words, we will manually specify a k-point list.\n", "\n", "Below we write the inputs, prepare the folders and run the SCF and NSCF calculations." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " \n", " \n", " -*#*- ...............- \n", " .+*= .+%*-=%%: .=#*- -===============-:.\n", " :*%=*%%- *%* #%* :+%+-%%+ .:. -=. :==-. \n", " -%S -%%*: :#%. -%%-. -##: #%* -=. :==- \n", " .. .%S: +%%%%*. :*%%%#= %%= -=. :==- \n", " :=#%%*- .#S- .. .%%= :*#*: -=. :==- \n", " -%S:.=#%%*==%# *%%=::=##-+%%. . .=-. :==- .= \n", " :%%- .-+++: -+##*=. =%S :-::==: .==- --.\n", " *%# #%+ -=--:. .----:. \n", " :%%- -%S. \n", " .-=*####SS%#########: -###########*+: +#######= =%SS####+::. \n", " =+#**%SSSSSSSSSSSSSSSS= =SSSSSSSSSSSSSSS= #SSSSSSS*. +SSSSSSSS%%%%- \n", " *%% =SSS.. .SSS= SSS=. .:+SSS+ -SSS:. .-::. .SSS+. #%* \n", " #%#. =SSS. *S# *S%: SSS= %SS%. .SSS= #SSS%. :SSS: =%#. \n", " *%%: =SSS#*#SSS- SSS= .+SSS+. %SS*.=SSSSS+ +SS%.:+%+ \n", " +%%:=SSSSSSSSS- SSSSSSSSSSSSS=. +SSS:SSS%SSS:%SS*-#%= \n", " ....#S==SSS:..SSS: =+- SSS%######+=. -SSS%SS%.#SS%SSS==%%=. \n", " .:+##%%#*- =SSS. ::. :SSS. SSS=. SSSSSS:..SSSSSS: :*%%%*=- \n", " #%+. -+#SSS*+++++++*SSS: .=+SSS#++++: %SSSS+. =SSSS%. :-+#%%+ \n", " #%* .SSSSSSSSSSSSSSSSSS: *SSSSSSSSSSS. +SSS%. %SSS*. . .%%= \n", " =%S :::::::::::::::::. .:::::::::. :::. :::. =+=-.#%+ \n", " -%S: =+===#S: \n", " ==*------------------------------=========+++++++++++++++++++++++========++-+## \n", " =+++++++++++*******++++++++++++++++========------------------==========+++++++- \n", "-- -- -- -- -- -- -- -- -- -- -- -- Structure Info -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "1\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si.upf\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si_r.upf\n", "pseudo found at pseudodojo : ONCVPSP-PBE-FR-PDv0.4/Si_r.upf\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "import EPWpy\n", "from EPWpy import EPWpy\n", "from EPWpy.utilities import parse\n", "size=3 # Supercell size given by size X size X size\n", "T=0.0 \n", "filename3 =str(os.getcwd())+'/'+f\"{prefix}/zg/ZG-scf_{dim1}{dim2}{dim3}_{T}0K.in\"\n", "atomic_labels, atomic_positions, cell_param, add_param = parse.parse_data(filename3)\n", "\n", "si_zg=EPWpy.EPWpy({'prefix':prefix,\n", " 'restart_mode':'\\'from_scratch\\'','ibrav':2,'nat':2*dim1*dim2*dim3,'calculation':'\\'scf\\'',\n", " 'atomic_species':['Si'],'mass':[28.0855],\n", " 'atoms':atomic_labels,'atomic_pos':atomic_positions,'ntyp':1,'pseudo':['Si.upf'],\n", " 'ecutwfc':'20','ecutrho':'80',\n", " 'celldm(1)':10.262*dim1,'verbosity':'high','atomic_position_type':'angstrom',\n", " 'pseudo_auto':True },env='mpirun',system='si_zg')\n", "si_zg.run_serial = True\n" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: scf -- -- -- -- -- -- -- -- -- -- -- \n", "Running scf |████████████████████████████████████████| in 1:55.5 (0.01/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "si_zg.scf(system={'nbnd':'119'},electrons={'conv_thr':'1E-6','mixing_beta':'0.7'},\n", " kpoints={'kpoints':[[2,2,2]],'kpoints_type': 'automatic'},name='scf')\n", "si_zg.prepare(1,type_run='scf')\n", "si_zg.run(4)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " \n", " \n", " -*#*- ...............- \n", " .+*= .+%*-=%%: .=#*- -===============-:.\n", " :*%=*%%- *%* #%* :+%+-%%+ .:. -=. :==-. \n", " -%S -%%*: :#%. -%%-. -##: #%* -=. :==- \n", " .. .%S: +%%%%*. :*%%%#= %%= -=. :==- \n", " :=#%%*- .#S- .. .%%= :*#*: -=. :==- \n", " -%S:.=#%%*==%# *%%=::=##-+%%. . .=-. :==- .= \n", " :%%- .-+++: -+##*=. =%S :-::==: .==- --.\n", " *%# #%+ -=--:. .----:. \n", " :%%- -%S. \n", " .-=*####SS%#########: -###########*+: +#######= =%SS####+::. \n", " =+#**%SSSSSSSSSSSSSSSS= =SSSSSSSSSSSSSSS= #SSSSSSS*. +SSSSSSSS%%%%- \n", " *%% =SSS.. .SSS= SSS=. .:+SSS+ -SSS:. .-::. .SSS+. #%* \n", " #%#. =SSS. *S# *S%: SSS= %SS%. .SSS= #SSS%. :SSS: =%#. \n", " *%%: =SSS#*#SSS- SSS= .+SSS+. %SS*.=SSSSS+ +SS%.:+%+ \n", " +%%:=SSSSSSSSS- SSSSSSSSSSSSS=. +SSS:SSS%SSS:%SS*-#%= \n", " ....#S==SSS:..SSS: =+- SSS%######+=. -SSS%SS%.#SS%SSS==%%=. \n", " .:+##%%#*- =SSS. ::. :SSS. SSS=. SSSSSS:..SSSSSS: :*%%%*=- \n", " #%+. -+#SSS*+++++++*SSS: .=+SSS#++++: %SSSS+. =SSSS%. :-+#%%+ \n", " #%* .SSSSSSSSSSSSSSSSSS: *SSSSSSSSSSS. +SSS%. %SSS*. . .%%= \n", " =%S :::::::::::::::::. .:::::::::. :::. :::. =+=-.#%+ \n", " -%S: =+===#S: \n", " ==*------------------------------=========+++++++++++++++++++++++========++-+## \n", " =+++++++++++*******++++++++++++++++========------------------==========+++++++- \n", "-- -- -- -- -- -- -- -- -- -- -- -- Structure Info -- -- -- -- -- -- -- -- -- -- -- -- -- \n", "1\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si.upf\n", "https://raw.githubusercontent.com/PseudoDojo/ONCVPSP-PBE-FR-PDv0.4/master/Si/Si_r.upf\n", "pseudo found at pseudodojo : ONCVPSP-PBE-FR-PDv0.4/Si_r.upf\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "import EPWpy\n", "from EPWpy import EPWpy\n", "from EPWpy.utilities import parse\n", "\n", "filename4 =str(os.getcwd())+'/'+f\"{prefix}/zg/equil-scf_{dim1}{dim2}{dim3}.in\"\n", "atomic_labels, atomic_positions, cell_param, add_param = parse.parse_data(filename4)\n", "\n", "silicon_equil=EPWpy.EPWpy({'prefix':prefix,'restart_mode':'\\'from_scratch\\'','ibrav':2,'nat':2*size*size*size,\n", " 'calculation':'\\'scf\\'','atomic_species':['Si'],'mass':[28.0855],\n", " 'atoms':atomic_labels,'atomic_pos':atomic_positions,'ntyp':1,'pseudo':['Si.upf'],'ecutwfc':'20',\n", " 'ecutrho':'80','celldm(1)':10.262*dim1,'verbosity':'high','atomic_position_type':'angstrom',\n", " 'pseudo_auto':True },env='mpirun',system='si_equil')\n", "silicon_equil.run_serial = True\n", "\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: scf -- -- -- -- -- -- -- -- -- -- -- \n", "Running scf |████████████████████████████████████████| in 28.6s (0.05/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Run SCF for the equilibrium structure\n", "silicon_equil.scf(system={'nbnd':'119'},electrons={'conv_thr':'1E-6','mixing_beta':'0.7'},\n", " kpoints={'kpoints':[[2,2,2]],'kpoints_type': 'automatic'},name='scf')\n", "silicon_equil.prepare(1,type_run='scf')\n", "silicon_equil.run(8)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Question}}$\n", "\n", "After you run both SCF calculations, the folders si_equil and si_zg will contain the SCF for both supercells. Compare the energies and band gaps. Are they the same?\n", "\n", "\n", "Next we run the NSCF calculations for both systems:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: nscf -- -- -- -- -- -- -- -- -- -- -- \n", "Running nscf |████████████████████████████████████████| in 1:34.9 (0.01/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Run NSCF for the ZG structure\n", "si_zg.nscf(system={'nbnd':'135','nosym':'.true.','ecutwfc':'20'},electrons={'conv_thr':'1E-6','mixing_beta':'0.7'},\n", " kpoints={'grid':[2,2,2],'kpoints_type': 'crystal'},name='nscf')\n", "si_zg.prepare(1,type_run='nscf')\n", "si_zg.run(8,type_run='nscf')" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: nscf -- -- -- -- -- -- -- -- -- -- -- \n", "Running nscf |████████████████████████████████████████| in 2:15.9 (0.01/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Run NSCF for the equilibrium structure\n", "silicon_equil.nscf(system={'nbnd':'135','nosym':'.true.'},electrons={'conv_thr':'1E-6','mixing_beta':'0.7'},\n", " kpoints={'grid':[2,2,2],'kpoints_type': 'crystal'},name='nscf')\n", "silicon_equil.prepare(1,type_run='nscf')\n", "silicon_equil.run(4,type_run='nscf')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, we prepare folders, create inputs and run the calculations for the real and imaginary part of the dielectric function with epsilon_Gaus.x from Quantum Espresso:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "obtaining nscf and ph attributes\n", "-- -- -- -- -- -- -- -- -- -- -- Calculation: eps -- -- -- -- -- -- -- -- -- -- -- \n", "on 1: running: mpirun -np 1 /mnt/storage/sabya/For_video/epwpy/build/q-e-EPW-5.9s/bin/epsilon_Gaus.x -in eps.in > eps.out\n", "Running eps |████████████████████████████████████████| in 2.8s (0.70/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Run dielectric function calculation for the ZG structure\n", "si_zg.reset()\n", "si_zg.verbosity = 2\n", "si_zg.eps(inputpp={'calculation':'eps'},energy_grid={'intersmear':'0.03','wmin':'0.2','wmax':'4.5','nw':'600'},name='eps')\n", "si_zg.prepare(1,type_run='eps')\n", "si_zg.run(1,type_run='eps')" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "-- -- -- -- -- -- -- -- -- -- -- Calculation: eps -- -- -- -- -- -- -- -- -- -- -- \n", "Running eps |████████████████████████████████████████| in 1.8s (1.57/s) \n", "\n", "-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- \n" ] } ], "source": [ "## Run dielectric function calculation for the equilibrium structure\n", "silicon_equil.eps(inputpp={'calculation':'eps'},energy_grid={'intersmear':'0.03','wmin':'0.2','wmax':'4.5','nw':'600'},name='eps')\n", "silicon_equil.prepare(1,type_run='eps')\n", "silicon_equil.run(2,type_run='eps')" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAEMCAYAAAAF2YvKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAABR60lEQVR4nO2deXxcddX/32e27FuTplu60wKltKVUkE0oi4BsAo/IJosoivAA6uMjKIsPioAbiKKIiiiy/ihrqStQQfZSW2gpXeiabkmzZzL7/f7+uDPJJJkks9y5M0m+79crr0lm7tx75mbmfuac8z3niFIKjUaj0WjSwZFrAzQajUYzfNEiotFoNJq00SKi0Wg0mrTRIqLRaDSatNEiotFoNJq00SKi0Wg0mrTRIqLRaDSatNEiotFoNJq0ceXagMEQkTnA94Am4CWl1FNDPaempkZNmzYty5ZpNBrNyOK9997bp5Qam+rzbBcREXkQOB1oUErNjbv/FODngBP4nVLqTuBU4BdKqddE5HlgSBGZNm0aK1asyI7xGo1GM0IRkW3pPC8X4ayHgFPi7xARJ3AfpmjMAS6IeiEPA+eLyI+Bapvt1Gg0Gs0Q2C4iSqlXgeY+dx8GbFJKbVZKBYHHgbOUUg1KqauBG4B9Npuq0Wg0miHIl5zIJGBH3N/1wOEiMg34DlAC/HigJ4vIlcCVAFOmTMmelRqNRqPpRb6ISEKUUluJisMQ2z0APACwaNGifm2JQ6EQ9fX1+P1+y23UjDwKCwupq6vD7Xbn2hSNJu/JFxHZCUyO+7suep8l1NfXU1ZWxrRp0xARq3arGYEopWhqaqK+vp7p06fn2hyNJu/JlzqRd4FZIjJdRDzA+cDzVu3c7/dTXV2tBUQzJCJCdXW19lo1miSxXURE5DHgTWB/EakXkSuUUmHgGuBvwDrgSaXUWouPa+XuNCMY/V7RaJInF6uzLlBKTVBKuZVSdUqp30fvX6aUmq2UmqmUut1uu7KN0+lkwYIF3T933nln2vs68sgjAdi6dStz55qlNitWrODaa68F4Hvf+x4/+clPUtqXHdxzzz10dXVZtr8f/vCHlu1Lo0mHQNDg8X+0s3rj6PVc8yUnMuIpKipi1apVluzrjTfe6HffokWLWLRoUdL7CIfDuFyuhPvKFvfccw8XX3wxxcXF/R6LRCI4nc6U9vfDH/6Q73znOyk9J/a6NZpU+dOyNtZ8HOC7l1dTUWq+V5/5VycPPNMKwJc/W8n5J5WNOk82X3Iio5a//vWvHHDAASxcuJBrr72W008/HejvTcydO5etW7cCUFpa2m8/y5cv734uwOrVqzniiCOYNWsWv/3tb7u3OeaYYzjzzDOZM2dOr331ff4111zDQw89BJhdAG688UYWLFjAokWLWLlyJSeffDIzZ87k/vvv72eL1+vltNNOY/78+cydO5cnnniCe++9l127drF48WIWL17cfexvfvObzJ8/nzfffJNp06axb59ZDrRixQqOO+44ADo7O7n88ss5+OCDmTdvHkuWLOGGG27A5/OxYMECLrrool5eGcBPfvITvve97wFw3HHHcf3117No0SJ+/vOf895773Hsscdy6KGHcvLJJ7N79+7k/lmaUcvW3SEeWtrGinV+fhMVDYC/v+Vl9hQPiw8t5rfPtnL7H5rwBYzcGZoD9Fcym4hd8GLceOONnHXWWXz5y1/m5ZdfZr/99uPzn/+8Zcd7//33eeutt/B6vRxyyCGcdtppAKxcuZI1a9akvPJoypQprFq1iq9//etcdtllvP766/j9fubOnctXv/rVXtv+9a9/ZeLEibz44osAtLW1UVFRwc9+9jNeeeUVampqAFNsDj/8cH76058Oeuzvf//7VFRU8MEHHwDQ0tLCueeeyy9/+ctu7y4msAMRDAZZsWIFoVCIY489lueee46xY8fyxBNP8N3vfpcHH3wwpfOhGV2s/MiPQ8JcecI/qXG+iH9jCx5HM/ecUY/XOJjxUz/LvKknce8z8P6mAF88o4ITDyvB5Rz5XsmoE5Ff/r8WPq4PWrrPmXUervlc1aDbJApnrVq1iunTpzNr1iwALr74Yh544AFLbDrrrLMoKiqiqKiIxYsX884771BZWclhhx2W1tLVM888E4CDDz6Yzs5OysrKKCsro6CggNbWViorK7u3Pfjgg/nmN7/Jt7/9bU4//XSOOeaYhPt0Op2ce+65Qx77n//8J48//nj331VVg5/rRMQEev369axZs4aTTjoJMMNoEyZMSHl/mtHBc//q4G9ve2lrb+Wez32JueNfZ593HA0tM6ismMGrm+Zx3Jw1SOP/cFYdfPrahfxt7ck8+sKneWjpLM49voxTjyyltMiBUmpEhrpGnYgMF1wuF4bR4xanuuS075s19ndJSUlaxysoKADA4XB0/x77OxwO99p29uzZrFy5kmXLlnHTTTdxwgkncMstt/Q7ZmFhYa88SLwNqb7eoeyPvW6lFAcddBBvvvlmSvsfLbR0RCgpdOBxj7yLXap4fQa/WtLCuNLN3HzydcyuXQvjH+CJf57Ds8u9nH50Kc+92smRJ02Coq3Q8TRFHUv47Jw7+OycO9jVfiD/WPdp7vzVQpr8C9iwq5KZdW4WHVjEIbMLmDuzgKKC4Z9RGHUiMpTHYCcHHHAAW7du5eOPP2bmzJk89thj3Y9NmzaNpUuXAmYIasuWLSnt+7nnnuPGG2/E6/WyfPly7rzzTjZs2DDg9lOnTuXDDz8kEAjg8/l46aWXOProo9N6Xbt27WLMmDFcfPHFVFZW8rvf/Q6AsrIyOjo6usNZfZk2bRrvvfcep556KkuWLOm+/6STTuK+++7jnnvuAcxwVlVVFW63m1AohNvtZty4cTQ0NNDU1ERpaSlLly7llFNO6XeM/fffn8bGRt58802OOOIIQqEQGzZs4KCDDkrrtY4ktuwK8uXb9zBxrIvf3DCeosLhf4FLhkhE8eslLZSVOLn0tIru+99evZlrj72JU+csIaKKCNY+TWHlmZy7OMwzr3h57tVO9qtzM6bCCcyE6m+ZP6Ed0PE0E4uWcEn5vQhmE43O0AS2txzI+t1TeevVqTyzdBquov2oHDOTGXUl7FfnYfpENyVFw+u8jzoRyRV9cyKnnHIKd955Jw888ACnnXYaxcXFHHPMMXR0dABw7rnn8qc//YmDDjqIww8/nNmzZ6d0vHnz5rF48WL27dvHzTffzMSJEwcVkcmTJ3Peeecxd+5cpk+fziGHHJLW6wT44IMP+Na3voXD4cDtdvPrX/8agCuvvJJTTjmFiRMn8sorr/R73q233soVV1zBzTff3J1UB7jpppu4+uqrmTt3Lk6nk1tvvZVzzjmHK6+8knnz5rFw4UIeeeQRbrnlFg477DAmTZrEAQcckNA2j8fDU089xbXXXktbWxvhcJjrr79+VItIY2uYzi6DZW94MRTUN4R54p/tXHZ6Za5Ns4UVH/l5enknAJ/9VAkVRY2EO17msKL/xnNAF1L1Fdw1N+N2jQNgfLWLS0+v4M9/aeO8E8v779A9GcZcB2OuQyLtEPgP+FdS6n+POaVrOHDcO4jq7N48YjjZ1TaFHetnsvStGbQEZoBrMp7iKZSVTWFcbTUTa9yMrXJSXuLIu5CYKNWv1dSwZtGiRarvPJF169Zx4IEH5sii5Fm+fDk/+clPuj0QTe4YLu+ZTFFK8YXv7WZXoxmS/NQhRUQMWL3Bz2M/mNTrW/Hrq7v4/fNtfPGMCo5e0H+Z9nDljy+28ccX25hStZG7L7iZKpcZ6ly/92D2Fv6BTx1+aMLnBUMqvbCfUhBpgOAmCG1CBTYQ9H5ExL+eAjbhlECvzX3BYho7J7DPO47mrgn4IhMISx24J+IpnExx6UQqK2upHVPC2CpX2qFIEXlPKZV8nUAU7YloNKOYLbtC3QICcNpRpZQVO3h9tY8XX+/s9U37oRfb2Lo7xO1/aOJ333UzqXYENKgMbmGS/Jr/O2Mbn5zyLGFVjFFzB3c8NpEu54nc/rWJAz417byRCLjGmT8chQDdWUYVgfBOMyQW3kHQtwN/xw6KnPVMLd3JLN6i2LUXpyPcb7dde0po2lJFV2gMSopwOj243B6crjIKiiZSXj4WX6icpW+62NZYSnHxGCrKx1A9poZxA4SYk0GLSB5x3HHH9QrjaDTZ5s0PfAA8+cOJKGBspXlJOGiGhxdf7+RzJ5jFcx1dBpt3hjjj6FJefs/Ljx5u5u6v1+JwpHchfeP9Lh5e1s7/XjKG6RM9Vr2c5Im0QuNNqNYHOH56mECkktV7z+HhFd/huov256V1e/j2JQlCVdlGnOCeYv4AnnLwjOuzjTJMTyZUTyhQj7dzN12djQSNBiKOfeBqgkgXwaAfv7+DEs9WynkZt9GOG7hwnrUmaxHRaEYxb37gY/+pHmoqe18KPnNUKQ8/t4rtG19j6vQjWb+tBqXg2IXFzJnh4a4/NfP08g7+6/j0LrR/fLGNjTtC3PtEC3d/ve9VMssYAdjxafCvxFf4RS6770tcdPqBtAcMPtjaxuurTWGdt1+hvXYlizjANR5c43EXLaKyEioH2DRiKBpbIqyuD7Jpk49wqI3TPhlifFWnKaSRFrp8zbS3NWGObkodLSIazSilsSXMh1uCfPGMnhVJKAX+tzmp7hd8+pLHcRoGbC7E1/wocCjTJ7k5ZP8C/rWyi98/18bhBxUxeVxqYS2vz2DjjhAlRcLqjQG27wkxZbyNobGmO8D/LkxawmtrP80+bzMHTi+g02cuEf/D0jbqal2Mr06tDU8+4nQI46tdjK92cdT8YhJNGS8uh+JxkK6IDK+1ZBqNxjL+Hf3GfezCaJLc+wpsPQS2HYHL9zxv77qKm5YtQbln88mqizhqvzep7PoGsut8vnk+uF3wo4ebiBipLc6pbwgB8MUzKgF4bZV1TTmHJLQdmu+CsvOg7Bz+tbKL2jFOZk12c+C0nrDaMQuK824VVL6iRUSjGaWs2uBnfLXT9CRafg07jgejA8Y/ADN30VV+F29sWsQGXqC5ayLf/8z5SMu90PEk1YFvc815Y1i7OcjTr3SkdNz6BjMpfMj+hRw4zcO/V/my8fL6oxTsvc78vfZHdHYZrFjn59hDTMEoKnBw65dqWLyomAtPzkE+ZJiiRcQGnnnmmV5t4BcsWIDD4eAvf/kLABs3buT0009n5syZHHrooSxevJhXX301qX23trbyq1/9yjJbt27dyqOPPmrZ/jT5iVKKDzYFOHhmAbQ+CHu/BiWnw/QPoPLL4Czj8LlFOBzw3OslXPvkI2xovxwmPAJjvgVtf+Skees44uAifv98G3ua+q8WGojte0I4BCbWuDh6QTHrtwdpaE7++WnT8Th0Pgs1t4F7Km+830U4EueJYf5+8xdrhl3BXy7RZ8oGzj77bFatWtX987WvfY1jjjmGk08+Gb/fz2mnncaVV17Jxx9/zHvvvccvfvELNm/enNS+BxORvu1IkiFdEYlEIik/R5M79jZHGFO4hgvnfxP2fAlKToZJT4Gj54JaVuxg/qwC/vqml8bOCbjr7oeKC6H6u+CsQRq+wXXnVaKU4tG/tSd97HVbg0yb6MbjFg6bYyav/7Mhy/M4wo2w5xoo/CSM+QYAy6OhrPgwliZ1tIjYzIYNG7jtttt4+OGHcTgcPPLIIxxxxBHdDQ7BbPt+2WWX9Xvu2rVrOeyww1iwYAHz5s1j48aN3HDDDXz88ccsWLCAb33rW/3avQ/WIn3Tpk2ceOKJzJ8/n4ULF/Lxxx9zww038Nprr7FgwQLuvvtuHnroIa655pru559++uksX74c6N/K/c9//nO3fV/5yle0sOQbkVYzqbz7Sor2LuZ3F57K5KL/B5VfhUlPg6Og31OOiRYVnviJ4p6luM4K89u87zVqPS/wmaNK+csbnUl5E4ahWLclwEEzzGNNn+imvMTB6o2BIZ6ZIY3fBaMdJvwexElze4R3P/Rz3EKd+8gULSI2EgqFuPDCC/npT3/KlCnmOvC1a9eycOHCpJ5///33c91117Fq1SpWrFhBXV0dd955JzNnzmTVqlX8+Mc/BsxeWz//+c8HbXMCcNFFF3H11VezevVq3njjDSZMmMCdd97JMccc0932fTBirdxXr15NdXU1TzzxBK+//jqrVq3C6XTyyCOPJPW6NDYQboRth0Pjd8yQTngP9//7u4Sm7oDxv+rlgcRzxtGl3H19Ld++tM+qnsovgWc2NN3F508sJ2LAX9/yDmnGtj0hvH7VLSIOhzBnuof126ztrN0L37vQ9juzFUmBOUfnqZc7iBjmUmZNZoy+Jb57rwf/Kmv3WbgAxt0z5GY333wzBx100KBzQ84++2w2btzI7Nmzefrpp3s9dsQRR3D77bdTX1/POeec091Cvi/JtHvv6Ohg586dnH322eZLKEx9TXx8K/eXXnqJ9957j0984hOA2SustrY25X1qskB4H+y+GEJbYcpyKD6Wu37VwN7mCF8t6r/kMx6nU5g/O8F7Q1xQeTU0XMf4CRtYMKua5e91cclnKvpvG8fazaZYHDSjJ4Q0s87DOx+2p99GZDBUGPZ8xayrqDY7Sdc3hFjycjsnf7KEKSkuT9b0R3siNrF8+XKWLFnCL3/5y173H3TQQaxcubL772eeeYaHHnqI5ubmfvu48MILef755ykqKuIzn/kML7/8csJjxbd7z7Sl/GDPj2/lrpTi0ksv7c77rF+/vjtspskh/lWweTZ4/w7j7oPiYwHYuCPErMkZ5gLKorNgOp/j8LlFbN0dorFl8JDW2s0BKksdTKzp+f66X50HwzCnB1pOw/+YDRBrfw5Oc8XVklc6EBG+dFal9ccbhYw+TyQJj8FqWlpauPzyy3n00UcpKyvr9diFF17IHXfcwfPPP9+dF+nqSrxufvPmzcyYMYNrr72W7du38/777zN//vzuzr+JGKhFellZGXV1dTz77LN89rOfJRAIEIlEutu1x5g2bRq/+tWvMAyDnTt38s477yQ8zgknnMBZZ53F17/+dWpra2lubqajo4OpU6emero0VtJ4g1nhPO0/pscMtHVGaGqLMLMuw2/h7klQuAg6nmPhAd8E4IOPAxy/aODLytrNZj4kPg8xISooe5vDzJ5iYZK7+W5o+TlUXQflnwPMnMxr/+nik3MLqa4Y/sWE+YD2RGzg/vvvp6GhgauuuqrXMt8nnniCoqIili5dyv3338+MGTM44ogj+MEPfsBNN93Ubz9PPvkkc+fOZcGCBaxZs4ZLLrmE6upqjjrqKObOncu3vvWtfs9xu93dLdJPOumkXi3SH374Ye69917mzZvHkUceyZ49e5g3bx5Op5P58+dz9913c9RRRzF9+nTmzJnDtddeO2D+Zs6cOfzgBz/g05/+NPPmzeOkk07Ss8tzje9t8P4Nxvxvt4AA7N5negsTx1rwHbL0LPC/zdSafTgcg3sTbZ0R6hvCzJnRO4Efu5g3t1m4EKPjOWj4hukt1faMX964I0hzu8FR80ZOF+Kco5QaUT+HHnqo6suHH37Y7z6NZjCG/XvGMJTaepRSG8YqFeno9dDL73aqxVdtU5t3BjI/ju99pdahVMtv1CXf26luur9hwE1XrfepxVdtU++s7ep1f9i/Tf39mbPVujcvVSrcnrlNgS1Kra9UassipSK9j/WnZa3q+K9tUy3t4cyPM8IAVqg0rrmjL5yl0YxUIm3Q8TR4ZkJwM/heh/G/A0fvFUgxT2R8tQUf/4K54J4OHc8xedzZhLs2wo7LzC604+4FcUNgLQTX09g4Gyjv6ZMV3gftf8K57zaOn92J0xGBfVUw7u707VHKLJwkAhOfAEdRr4ffWetj9mQPlWU6lGUVWkQ0mpGA4YXtR0NgTc99RZ+Cisv7bbqrKUxVucOa+d4iZkir5R4umHsV+1c9D95oSCu8E1QIvH8F4MRqmHf5JMb6T4T6FuhcBoSh+Hhuffb7nHHQ/RzOL6HySihIcyBY51Lw/gVqfwaeGb0e8voM1m0JcuEpuqWJleiciEYzEmj9vSkgEx+D8b+Fmv+DuufMpHofdjWGmWCFFxKj4guAmznVT/PPDZ+Hmduh+jvQ+QL4V8LYH8LUt3hhw+1sb1uAeJea94+5HqathikvEZT9WPLB/5pe056vmjMzUsXwQ8P14JkDVdf0e3hTfRBDwdwZ/YsqNekzajwRpZSuTNUkhRpuI6OVgtZfQeHhUH7+kJvv3hc2e2ZZReFCmLmVh5e18Pi/yjnl7Mkw9nazEt5Z0x1SemHtFKrLv8SiE/vXD7ldQmNHtZkE33MFtP4aqq5OzY7mn0JoM0z+pxlG68PH9aaHNLNOtzmxklEhIoWFhTQ1NVFdXa2FRDMoSimamprSKr7MGV3LIbgeJjw85KahsDmkaEKNxR9990QMZzH+YDuGocyJh+7JvTZpao2w/wBLeN0uCEUww28dS8yi4OBGcM+AwvlmaG6wz25wKzT90FyNVXJCwk027wxSWepgTLkOwFjJqBCRuro66uvraWxszLUpmmFAYWEhdXV1uTYjeTr+H0hxT/HfIOxrjWAoi5LqfSiM5lgCQUVRYe8LfjiiaOkw+k1QjOFyCuGwMoVi4qOw+3Jo+QUQDWuVftZsECkJEuLKML0XcZq5kAGobwgzeZxbf5G0mFEhIm63e8g2IBrNsEQZZi+s0lP7rURKRHO7WYtRVW796qQij3lx9gUVRX0cuaZoDchABX5ulxCKRMOIzgqoexpUBCLN0Ho/7LsF2p8wuwj3pfXX0PWyOQclOps8ETsbwyw6cBh5mMME7ddpNMMZ/9sQ3g2lZye1eUtURMZkQ0QKoiIS6J8U39dqHremMvFxXU5MTyQecYJrrNl63jMbWu/r/8T2JeagqZJToOJLA9rmDxo0tUWsKbDU9EKLiEYznOl4GnBD6WlJbd7SEfNErP/ox8JZ/kD/hQkxERk7kIi4hPBABevigPJLwPcGhOp77u/8K+y6AIoOg4lPDpoz2dtk7nyi1bkgjRYRjWbYohR0PAUlJ4GzMqmntHSYXkJVFortCqPhLH8wkYiYBY4DeSIelxDq64nEE8v3dDxj3oYbYNfnoeAgqFsGzrKBn0tPOK1G98uynLwWERE5RkTuF5HficgbubZHo8krAivN9u5ln0v6Kc1tEcpLHLic1ieXC4cIZ7ldUF6S+JLjcgrhyCAiUnCAWf/RsSS6w++ZBZYTH0tKQGMiMkaLiOXYLiIi8qCINIjImj73nyIi60Vkk4jcAKCUek0p9VVgKfBHu23VaPKa1ocAN5SdOdSW3bR0RKgqy87HvmiwcFZbhJoK54Aro9wuCIWHqNEpOxd8r4F3ObT+BiqvMsUlCZraB0/sa9InF57IQ8Ap8XeIiBO4DzgVmANcICJz4ja5EEh98LdGM1IJ7YC2B6H8AnCOSfppLR1GVlZmARREw1mBUH8haGqLUD3A8l6g2zOKDFaoXvZfgAE7FoOjHGpuTdq25rYIhQVCcWFeB1+GJbafUaXUq0DfiUuHAZuUUpuVUkHgceAsABGZArQppQYemqHRjCYMn5kPQEHN91J6akt7JGsi4og6GUYCIWj3GlQMEMoCM7EODJ4XKZxneh+OCpjwB3DVJG1bU3uE6iy97tFOvsjyJGBH3N/10fsArgD+MNiTReRKEVkhIit0QaFmxLP3GvC9ZVaoe1Krf2puj2QlqQ4QC1QlkoF2r0HZICLijpo0qIiAOQ9+diuUfTYl25rbIjofkiXyRUQGRCl1q1Jq0KS6UuoBpdQipdSisWPH2mWaRmM/nX8zw1jVN0D50BXq8fiDBr6AYkyWciISdUUS5TU6vMaASXXo8UTCg0/XTZumNu2JZIt8EZGdQHyjnbrofRqNJkakA/ZcCZ4DoPqWlJ/e0h5d3puli2m3J9JHQwJBg0BIDSoi7mhOJDTYCq0MaG7Xnki2yJfKm3eBWSIyHVM8zsdMpms0mhiN/wPhHTDl3+BIvX1Ha6e5QqmiNEueSFRF+spAu9cUr/KSgS/iPZ6I9SLiCxh0+ZVemZUlcrHE9zHgTWB/EakXkSuUUmHgGuBvwDrgSaXUWrtt02jylrY/QesDMOZ/oPjItHbR2TX0xdwS+uhAj4gkkROxcMx6jNjs9mrdvTcr2O6JKKUuGOD+ZcAym83RaPKf9qdg92VQfALU3Jb2bjqiIlJanGVPJA0RyaYn0l2lr3MiWUFLs0aTz/hXwu4vQNERUPd8WmGsGDERKc+SiDiGDGcN4om4spcTaffGwnhaRLKBFhGNJl9RCvZ8DRyVMOlZcBRntLsOb3Y9EQb0RMyL+KBLfLPoicTEsyxbr3uUky+JdY1G05f2x8xW7+P/YLZEz5COLoMCt+BxZ2coU8/qrN5CEBOvisES61nMicQ8ocFETJM++qxqNPmI4YfGG6HgEKi4xJJddnYNXvCXKbG+WH09kTavQaFncPHKtificEBJYXbEc7SjRUSjySWRJnOCXzxKQdP3Ibwdan9iztOwgPYuI6shnYEq1tu9Qx831jtryIr1NOjwGpQWOfRY3CyhRUSjyQWBD2H7p2FjDWyeBV3/js4HeR62HQ5NP4TyS6HkeMsO2dllZC8fQpzWJVidVT5EbUrMSwkmaN6YKe1dg1fLazJD50Q0GrtQIWj4NrQ/CpG9ZiPB6huh/UnYfqw5Hzy0FdwzYPxvoeIySw/f0WUwoTp7H/nY93yjjw50JHER7x5olQUR6UjCE9KkjxYRjcYumu6ElrvNluZFR5tt3F21MOYGaPwWBNZDzf9B+YUg1n80O7wGsyZnMydi3vbNiXR0GUwZN/jr6fZEEkxFzJSOLiNrVfoaLSIajT2oELT8HErPgkn/r/djznIY/5usm9Dhy01Yp8tnUFKUpCeSJRGZXKsvddlCy7NGYwfef5hJ9IorcnL4cEThD6js5kS6PZHeQuD1Dy0iBVnMiXQM0YZekxn6zGo0dtD+KDiqoPTknBzejoI7R2yJb9x9hqHo8qshl9c6nYLLab0nEjEUnT6dE8km+sxqNNnG8ELHs1D+ORBPTkyIFfxl9WKaICfSFZ23PpQnAqY3kmi0biZ4fQZK6Wr1bKLPrEaTbTqfB+U1E+Y5wg5PJNE8kS6fedxkZpsXeIRAcLAh66nTLZ7Z7lw8itEiohmehHZA/bmweQ403WUmrvOVtkfBNQmKjsmZCd7oxTwZjyBdEq3O8vqTP242PJH2LDed1GgR0QxHvC/D1kOg6+/gHAONN8D2E8zEdb4RWA/eZVD+Bcsqz9Ohyx/zCLJXtZ2oYt3ri4azkjiux+OwPLGe9aaTGi0imjxBBUENEcpQCpp/ATs+Dc5amLoCpv4bJvwZ/O/A1iMguNEeewcj0gz158CGKthyADhKYMx1OTWpyx+7mGc/J0KankihWyxPrHd1H1+3PMkWWkQ09qEU+FaY387j79t3O2wog80zwf/ewM9v/gk0XAulp8PUt6Fgf/P+iotg8ktgtMC2o8C/uv9zDb/pwYQbrX1NfVER2HkueF+Ess9B1ddhynJwjc/ucYcgdjEvzmI4KzZPxIiLZ6WSE/F4xPJiQ180sV9UoC912UKfWY09KAV7r4JtnzC/ne+5BiIdsOfLsO8mKD7J9ER2nArBrf2f3/E8NH4bys6DSU+Ds6z348VHmbPHpQC2HwfN90DX6xBYBy33web9YMcJsHk2+Fdl73U23QVdy2Hcb2DCAzDuZ1C4MHvHS5KYJ1JUkMVwVoIh696YB5SEJ1DoEcvbntgRxhvtaBHR2EPns9D6G6j8GlRdB633wcYKaPs9VN8EdS/A5L+bYa2dZ0C4oee5/vdh14VQuAgmPDRwbqFgfzO85ZkNDV+H7UfDljmw9xpwT4MJfzIHO+3+Yv/eHFbg/Qfsu9lsZ1JxqfX7z4Auv0FhgeB0ZP9i2jsnEg0nJeOJuK33RGJLjJPxhDTpoXsBaLKPUrDvFvAcCON+bvaFKj0LOp6GksVQdo65XcH+MOkpqD/NvPhXXmkKx97rwFkZne5XNPix3FNh6ltmI8PgejM/4Z5mjpcVMUVqz5eg6yUoOdG61xjaAbsuMF/j+Ad6lirlCV6/kd18CAOvzhJJzgPKhifi8xu4XT2t5jXWo0VEk326XoHAmqgXEX3LlSw2f/pSciJMWwkN/wtNPwIi4BwPk/8K7onJHU8EPNPNn76UXwyNN0Hzj3tEJLQdWh8wBaji8tSbHyoFe64CwwdTnwZHaWrPtwGfX2U9pJMgmkWXz6C4UJKa5VHgFgKWJ9aV9kKyjBYRTfZpuQ8cY6Ds88ltX3AQTH7RDGkFPoSiRdZdmB0FMOZaaPyOGSZzjYWth0Nkj/l45wtmziUVIelcaibSa39qhtLyEK/foDjLyWVJuDpLJV2bko1iw66AkdU8kEbnRDTZJrgFOp+DyivAUZjac121UHKc9d/sK78CUgL7boWdnwejHaathtp7TRFp+EZq+2u+C9zToepaa+20kC6/ojjLy1wTzRPx+pIPo2Wj2NAXUFkXz9GO9kQ02aX5p4Ajvy6wzjFQdTU0/8j8e+KjUDjP/AltgpZ7Ta+p+Kih9+V7B3yvQ+09WZkBYhVdPoMJY7NrX6Iuvsl08I1R4HEQCptNE61aAODzGxTplVlZRUu0JnsEt5irryq+AO66XFvTm7E/hPEPwuR/mKupuu+/A5zjoOm23tsH1oD3n/1XdTX/DBzlUPHF7NucAV0BOxLr/S/WXb7kczHZaAffFdA5kWyjz64mO4QbYOeZIG6ouTXX1vRHnFB5ef8VWo5iqLrKXK4b2mXeF1gLWxbCjpNMUYwR2gYdT5mryPrWreQZXX5l2zfyvquzkhWvbMxZ9/l1TiTbaBHRWE9ol1nYF/zYXJbrnpJri1Kj7HxAQUd0AmHLL81YjWsK7LvNrEoHaP45IPkVqkuAUiql3EQmOKT36ixfIHnxysZ0wy6/0tXqWUafXY11BD4yL6zbDjNDWXVLoeT4XFuVOgX7Q8Eh0P6YWVfS/oQ5F33c3RDeAZ3LINICbb+D8vPAPTnXFg9KMKSIGDZVbUvvFmj+gEGhO8lwVlRErEyu+wKGrlbPMvmbCdQMLzqehp3nAREoOBjqXoTC+bm2Kn3Kz4PGG6HlN2ZPrvKLzdCXayK0/tps+Gh0wJhv59rSIfHa0XwxitDbE/EHFYVJegKxnIhVtSJKKTMnoj2RrKLPriZzVBD2fA0K5sPMreZy2eEsINBTRd9wrVnsWHKimd+p+BJ4/wJNP4Cyc80VXXmOL9o/qsgOEZGenEgobHpAsTDVUHR7IhaJSDCkMAz06qwso0VEkzmdyyCyF8Z+36z6zrOWH2nhmQ0lp5m/V11tCghA1VdNb8RZC2PvzJ19KdDjiWT//yJxOZFYbqMgVRGxKJyl+2bZgw5n5QvKgMb/NdtveA6Acb+AosNzbVVydD4Pjioo+XSuLbGWiX80W9MXx63gck2AGR+b3pezPHe2pUBPJ1s7PJEeV8QfrT5PNrHdE86ypmo99rr16qzsoiU6X2j9lVmYV3wihPeaczE6/5Jrq5LD+xIUL87rYru0cFabwti3a7CjcNgICMTNErHDE6GnYt0f9QSSD2eZ59mqcJavu/29vsxlk7w+uyJynIi8JiL3i8hxubYnaxhe2Pd980I8aQlMf99MTu+6MPFsjXwivBvC26E4d/PDNYPjt3EwU3xOJOVwltvacJYvoGeJ2IHtIiIiD4pIg4is6XP/KSKyXkQ2icgN0bsV0AkUAvV222oJ/tXQvmTwbVrug0gD1PzA/BQ6K8wmgITNQU7ZmH1hFf7/mLd5MHhJk5ie6X72LPGNEfMocpVY1zkRe8jF2X0IOCX+DhFxAvcBpwJzgAtEZA7wmlLqVODbwP/ZbGfmhPfA1gWw67+g/fHE20Q6zB5OJSdD8ZE993umQ83t4P0rdDxhi7lpERORggU5NUMzMLFv5Mkutc0Eoec7T+y4KedErPJEdE7EFmwXEaXUq0Bzn7sPAzYppTYrpYLA48BZSnWXLbUABTaaaQ3tcRf/1gcSb9P8I4g0Qc1t/R+ruhoKP2EOZYq0ZMfGTAn8B9wzh1WOYLThT9EjyARLVmdZ5YnonIgtpHV2RaQk6j1YxSRgR9zf9cAkETlHRH4DPAz8chB7rhSRFSKyorGx0UKzMsS7DDwHQfV3oetfptcRQxnQ+F2z3qD8C1B0WP/ni9OckhdpgsYb+j+eD/j/A4WH5NoKzSD4AwqXE9wum0QkmllPNZzlcgpOh5VLfHVOxA6SEhERcYjIhSLyoog0AB8Bu0XkQxH5sYjslw3jlFJPK6W+opT6vFJq+SDbPaCUWqSUWjR27NhsmJI6SoFvBRQdCUXHAIZZ5Ryj6Q5o+iFUXAHjfzvwfgoXwJjrTU+m699ZNjpFIm0Q2my2CNHkLb6AYYsXAr0r1tPxgAo81s0U8emciC0ke3ZfAWYCNwLjlVKTlVK1wNHAW8BdInJxBnbsBOIbENVF7xu+hLaA0WxO5Ss6HBDwvWk+Ft4drXg+zxQQxxCRupr/A9dU2HMlGJ1ZNz1pAqvMW+2J5DX+oH1NCEUkbnVW6rkYc7qhdTkRPV89+yS7sP9EpVSo751KqWZgCbBEJFbSmxbvArNEZDqmeJwPXJjB/nJPYK15WzAPnJXgmQO+N8z72v4Iyg9jb0+uuttRAhN+CztOgS0LoPx8M0TmKMqW9cnRvTJLi0g+4wsoCm1KLvfKiaRYJwKxOetWFRvqWSJ2kNQZjhcQEUnYcS6RyCRCRB4D3gT2F5F6EblCKRUGrgH+BqwDnlRKrU1mf3lLaKN565ll3hYdCb63zFxI+yPm354UooAlJ5lNDd2TzTDY3quttzlV/O+ZfaVc43NtiWYQ/AHD3uRyXJ2IwwGuFLKnBR6HpTkRvTIr+wzpiYjIk/F/AguAu9I9oFLqggHuXwYsS3e/eUdwo9kKxFlt/l10JLT91pxREVgD4+5LfZ+lp5g/e683a0tqfgDuiZaanRK+N8zXpclr/EFlW07E0afYsMgjCSceDoTpiVhXsa47+GafZM5wu1LqvOjP54B/ZtuoEUFwU48XAlB0hHm79xrAZeZD0qXySiAM3qWZWJgZoZ1mUj32ujR5i+3hrLicSKq1KVbmRLoCer66HSTzH769z9/fzYYhI47wzt7DijyzwVEJkX1mYaGrJv19ew4E1ySzZ5VdqAh0vGDWvhidPQWQpWfYZ4MmLWwNZ/XJiSRbIxKjwG3t6iydE8k+Q4azlFJbAESkRim1L5pM1wxFeFfv+d0iZs3Hvpug9seZ7VsEij9l35JfZcDOz0HnM+bfjlIw/FB0tDkFUJPX2BnOMivWe+pEUj1ugUdobrNudVZtlZXlbJpEpCLTD2bNipGG4QWjzZw7EU/552DGeig4MPNjFBxijmqNNGW+r6Fo+6MpIDW3w5RXzVGxpWfAhIezf2xNxvgCyrYEs4j0qhNJWUQs9ET0fHV7SKV3tw4uJkt4t3nbV0SsJLas1v+f3h6P1ShlFkYWLoLqG6NekO7YO5xIJzeRLhJXbegPpl7kWOCW7iLFTOnS89VtIZV3Vh63ks0z7BCRgoPN20CWV0IH/mMuV668emRMLBxlhCOKUBiKbAxnGXGrs9JLrGdeJ6KU0quzbCKVM6yvIMkS3mXeZlNEnLXmEuLguuwdA6BzKSBQelp2j6PJCt0Ff7lYnZVGYr2wwJo6kWBIYSg9X90OUhGRG7NmxUjDDhERMXMrgWyLyAtQ9Elw5UlPMk1K+FIcUZsxcdfsWJ1IKhR6hFAYIpHMhETPErGPpM+wUmrN0FtpAFNEpBAcFdk9jufA7Hoi4d3gXwElp2fvGJqs4kuj9UgmxM8TSScnEts+07yInq9uHynLdIY9skYH4V2mF5LtHILnQIg0Zm+FVme0gYCuBRm2+O2cakisYl11Hzs2Nz1ZrBIRPV/dPlI6wyLyW2C7iOwQkbdF5Lci8t9Zsm34EhORbBNbKpytkFbnUnBNgYK52dm/JuvYOdUQgGgX33BEETFSz8X0iEhmyXU9S8Q+Un1nfQqoU0pNBs4BngFKLLdquBPea09TQk9URLIR0jL84P0HlJ6uV2UNY2Lf6O2rEzGXcabTwRd6xC72/HSJTTXUOZHsk0qdCMDbQDXQoJTaidm2feQ0TbSKSFNP48Vs4p4KUgSBj6zbZ+tD0PYQOApBeaHsXOv2rbGddC/m6RLLicQS+ulUrEPm0w39AZ0TsYtUReQ3wL9E5PeYgvK+UqrNerOGMUpBpBmcY7J/LHGAZ3/rPJFQPez9Cqig+XfRp6B4sTX71uQEu8NZMU+kZzRumjkRizwRnRPJPqmKyJ+BB6LP+xowT0QKlVIzLbdsuGJ0AmFw2CAiYOZFYhMTM6X9CVNAZmwC5TObRupQ1rAm1TnnlqDSr0/ROZHhR6oiUq+UuiP+DhEZYrbrKCO2UsqOcBaYeZH2x8HoAkdxZvvq+qe5P4/+TjBSiIWFCtz2rc4yVHrz1c3tozkRi1Zn6ZxI9kn1DK8Skevi71BKBSy0Z/hjRJsc2xHOgmhyXUFwfWb7UcocMlV8rCVmafKDmCeSauV4ukh0ia8/05yIBXUier66PaQqIuOAr4rILhFZKiK3i8jnsmHYsMVuT8SqZb7h7WC0Q8H8zG3S5A2BkMLpsPFiKr1zGqnWicQS4b4McyJ6loh9pHqWfwR8ApgO3AJsAA6z2qhhTcRuT2QWiMdslJgJgQ/M21hjR82IIBBKvX9VJsRWZ8XCaOnmRDJdnaXnq9tHqjmRPwELoyGslSKyHbBhoMUwwm4REY/Zpt33Rmb7CUTDYQVzMrdJkzcEgjaLSLQBY8yTSLV3lsctiPQs0U0X3cHXPlI9y36llD/2h1JqH3CbtSYNc7rDWTaJCEDRkWaPK8M/9LYDEd4OjjJzhK9mxBAIGrYl1SGu2DCaE0lVwESEAk/mM0X0fHX7SFVENovIqX3u81hlzIgg0myOjxUbT0vRUebSXP/K9PcR2m62ONFLekcUZjjLvm/kMU8k3ToRgEILBlN1+XVOxC5SDWf9N/AXEfkC8BZwEPCx5VYNZyJN4LApqR6j6Ejz1vc6FB+Z3j7CO8A9xTqbNHlBIKjs9UQwRcQfUDgc4E71CoOZR8l0dZbPbzBujJ6vbgcpSbVSajdwKLAEGAu8D1yYBbuGL3ZVq8fjqgX3tMw9Efdky0zS5Ad250SILfENmfPVJQ3PtsDjyLjY0Jwrrz0RO0jqe4KIiIr2d1ZKRTBFZMlA24xqjByICEDBgvRXaBk+s6W8S3siI41ASFFSZGM4K3rrD6Q+SyRGoUU5EV2tbg/JvrteEZH/FpFeVxkR8YjI8SLyR+BS680bhtjVfLEvhQsguMGsXE+VcL15q8NZIw67w1mOaCv4dOarxyj0SEa9swwjOl9d50RsIdmzfAoQAR4Tkd0i8qGIbAE2AhcA9yilHsqSjcOLSCs4K+0/bnfl+sbUnxvabt5qERlx2F0nQjSx7g8qCtMUr0JPZjkRX8Ccr15WrEXEDpIKZ0WX9f5KRN4BVgM1gE8p1ZpF24YnRlv2x+ImwjPbvA1uhMIUq85jIuLSOZGRRiBkc2I9bp5IqoWGMQoLHPibwmnb0Okz8ymlNobxRjOpnuVLgEeBqTEBEZGfWW3UsEUFQfnBUW7/sT37mbfBDak/NxbOck2yzh5NXmB7sSExT8RI+7iFHsGXgSfS2RUVEe2J2EKqZ7kBs+3J0yKyPhrSsmGE3zDB6DBvcyEijlJwTYDQptSfG2k0iwwduiHzSCMnnki0TiSdGhHIPJylPRF7SXUV98XA/kqpgIhMBO4AMmzaNIKItJu3zhyEswBcU3tCU6kQ2QfOGuvt0eQUw1AE7e6dFT2UmVhP77gFGSbWu0VEeyK2kOpZ3oHZfBGl1C6l1KXAVyy3arhiRIc85sITAbPOI7wj9eeFG7WIjECC4fQ66WaCiGAolVFivbjQQSCkiETSExJvNJxVUqSX+NpBqp7IdcASEVmJ6YFMAryWWzVcMaKeSM5EZAp0vmDGE1Ip8ors04WGI5DuWSI5qVg30vZEYmEor9+gvCT1qvNOn+q1H012SbVi/UNgIfA4UATsAc7Kgl0AiMiBInK/iDwlIldl6ziWkWsRcU0xE/uRfak9T4ezRiR2D6QCoioSDWel6QGVRj2ImBikSkeXzonYScpnWSkVUEq9qJS6XSn1Y6XUzlSeLyIPikiDiKzpc/8p0WT9JhG5IXqsdUqprwLnAUelaqvtRKLhrJzlRCaat+HdyT9HKS0iIxS/zaNxwdSQsAHhSPpz3WMV9rFVVqnS2hGhvMSBU081tIVcSPVDmMWL3YiIE7gPOBWYA1wgInOij50JvAgss9fMNMi5J1Jr3kYakn+O8preixaREUcwB56IQ8AXnQWSaTgrliBPlZaOCFVl2guxC9vPtFLqVaC5z92HAZuUUpuVUkHMcNlZ0e2fV0qdClxkr6VpkGsRcY4zb8N7k39OOBr6co613h5NTgnkwBNBegZSpeuJxFZVpeuJtLQbVJXpDr52kUaj5qwwCXPlV4x64HAROQ44ByhgEE9ERK4ErgSYMiWHrTuMdsAJUpSb46fjicTyJ9oTGXHkIicicZ5IuqvCYuEsbwaeyKzJesyRXeSLiCREKbUcWJ7Edg8ADwAsWrQod52EI9GWJ7ka7OSoBNwQTkVEGs1blxaRkUYuPBGB7hqPtD2RTMNZ7TqcZSf5cqZ3AvFrTOui9w0vjHZw5iiUBaZ4uWohkkI4K6LDWSOVQJojajNBRHrmq6eZEykuNOespxPOCgQNvH5FVbkOZ9lFvojIu8AsEZkuIh7gfOD5HNuUOkZ77vIhMZy1KXoiOpw1Uun2RGwsNownXfFyOITyEgetnamLSENLBICxVXkdZBlR2P7uEpHHgDeB/UWkXkSuUEqFgWuAvwHrgCeVUmvtti1jjLbci4irNsWcSCPgzE3nYU1WiQ12SrdyPB0ccVeUdOtEACrLnLR2RFJ+XmNURGqrtCdiF7bLtVLqggHuX8ZwWMY7GJF2swliLnGOg8BHyW8faQFnVe7yOJqsEYx6Ih6bu/jGSDcnAlBV5qClI3VPZG+L2UK+doz2ROwiX8JZI4Nc50SgJyeS7KRioy2akNeMNPw5aHsSryKZiUhmnkhNhfZE7EKLiJXkQzjLWWsWDxqdyW2fq0mMmqwTDCocDnDZeD3t5YmkmViHqCfSnrqINDSHGVPuwGOncI5ytIhYidEBjrLc2uCKFhwmmxfJ1SRGTdaJzRIRG0OV0ssTySwn4vWr7pBcsjS0RHRS3Wa0iFiFikSnGuZYRGKrrGL1H0MRadUiMkIxB0PZ+408JlgOAXcG1/LYEt2WFENajS1hxo3RoSw70SJiFUa0I76U5NYOZ5V5G2lNbnujTYezRiiBkLI9rBPzRAoLMvOAKqPFgqnkRZRS7G2JUKs9EVvRImIVsRyEozS3djiiImK0JLe9DmeNCJrbIjzwTAsNzeHu+/xBw35PJHqbaYHjmLKYJ5L8Cq2OLgN/QDFWL++1FS0iVqHyRERiXkUkCRFRYVP8tCcy7Lnr4SYe/0cHt/2+Z5ZMMAeeSExFMsmHgJkTgdTCWbsaTQGdOFZ7InaiRcQqYuEsR47DWY4UwlndXYe1JzKcaW6PsGKdn8IC4cMtQXY2hIDMBkOlS0yyijL0RHrCWcl7IvUNpojUjXVndGxNamgRsYq8CWcVmF2EkwlnxYZoaREZ1ry3zo9S8O0vVAOwcr0fyI0nEqtYzzScVVTgoLBAUlrmW98QwiEwoUZ7InaiRcQq8kVEwAxPJeWJtPZsrxm2bNwRxOMWjp5fRHmJg/XbgoDpidg6GjeOTGpEYlSVOlJKrG/ZFWJ8jUvXiNiMFhGryJfVWWCGtLQnMmrYuCPIzElunE5h9hQP67ebIhIM2S8isRVZVlTJV5U7k06sK6X4YFOAuTMKMj6uJjW0iFhF3nkiSYiIkeOZ8BpL2Lo7xIxJZh5g2gQ39XvDGIYiEFT2tjwhLidSkPmlpbLMmXRi/f1NAVo7DebP1iJiN1pErELFEut5ICKOqtTCWbp31rCl02fQ1mkwqdYUkUm1LgIhxb7WiFmxbrsnYt5asbTYbMI4tIi0eyP8+M/N1FQ6OW5hccbH1aSGFhGr6PZE8iCc5UwxnKU9kWFLbFnrpOiy1slRMalvDOfGE4kezgrxqixz0tZhYBiDtz55/O/t7N4X5pYraizxgDSpoc+4VRidgORuvno8jsokw1mt0e21iAxXdjWay3ljIhKrkdjZECIQykXbE/O20IKLeVWZE0NBu3fgvIhhKP72tpej5xcxd6YOZeUCLSJWYXjNUFY+zOVwVpn5DjVEUjLSZi4EEL0kcrjS0D2Eyfwf1lQ6EempmbD7m7kRfctlWicCdM9JHyyktak+REu7wVHzdRgrV2gRsQqjMz9WZkG0f5bqKSYcCKNNh7LylA3bg9z1pyY2RFdaDURzewSPWygpMi/aLqcwptzJ9j2mh5LunPN06RnJa83qLBi84PDDLQEAFszSXkiu0CJiFaozP5Lq0JMoHyq5bnTmvuuwph+RiOL7v9/H397y8rNHm1GDDBhraY8wptzRq9lhTYWT7XtNT6S40N6PeCz0VFZszeosGMIT2RGkvMSh+2XlEC0iVhELZ+UDsU6+QyXXjTwSPk03H3wcYGdjmHn7FbBhe5Ade8MDbtvcbjCmvPcFtKbKye59sXCWvZ5IQ3Q87TgLxtNWJdH6ZFN9iFmTPbbOTNH0RouIVeTTBbnbE0lCRCRPbNZ08+YHPtwu+O/zzC8DsTYmiWhpj3SHfWKMrez5u8hmT6Sh2fQarJjpUVbswOFgwNYnSinqG0JMGa9zerlEi4hVGJ35sbwXeua8Gx2Db6fyyOZRwsqP/Kz5ODDoNms3BzhgWgEzJrmprnCybsvA2ze3R/p7IvEiYrMnEsOK6YIOh1BV5qSlM7GItHUadPkVE3WvrJyiRcQq8imcFctzDCUi+WTzKODdD338z70NXH/33l5zP+IJRxSb6kMcMNUM0cyY5GbL7tCA27Z5+4ezxlb2XFTtXp01bYJZp2JV/6rKMgct7YnDWTu7a2R0195cokXEKvIpNOSIeSJDrc7KoxDcKOD5V82CVMOAf77rTbjN1l0hgiHF/lM9AEyf6Gbb7hCRBAV3rR0RlIIx5b0/xvFJ5mKbPZF7/2ccj/9gomX7qypzDtiEMZb30V17c4sWEatQ3vwJDSXtiehwll0YhmLVBj+nH13KjEluVnyYOM/xUbQD7/5TTBGZOsFNKAx7mvp7Ls3Rb+iDhrNszomUFjmotSCpHqNykNYnsSR+rZ6pnlO0iFhFPn2rlyLAocNZecTW3SG8fsXcGR4WzCrgo23BhO081m8LUFbs6K48n1ht3u5p6n8hbY4mnPuKSPyM8VzlRKzC9EQSh7MaWyOUFTt0q5Mco8++FagIKH/+fKsXMb2RwUREBYFQ/oTgRjib6s28xgHTCphZ58EfVN3hmHjWbwsye0rPktXx0VBNom1jq5b6rs7yuIVLT6tg1mQ3LufwFxF/UOHz9xeSfS2RXivRNLlBi4gVGF3mbb5UrMPQIpJPDSMtpNNnJLzg5Jrte0I4HWZvq1jb9o939k6YB4IGm3eZSfUYNZVOXM7EIjKQJwJw6WkV/ObGCVa+hJwQG5Pb0tn/f9rYGukVutPkBi0iVqB85q0jD5ovxkhaREaOJ7Jua4D/umEnV/94L4FgfgnJ9j0hJo114XIKU8abIrJjb28R2VQfwjDoTqoDOB1C7RgXexKs5mpuj1BaJCN6kl/My0pUK9LYGtYikgdoEbGCbk8kj5rADSkieTT/xCIe/Ws7wZBi6+4QL7/XlWtzelHfEGbyOFM8igocVJU7+nkXsbG2+0/z9Lq/psLJvtZEORGjXyhrpFFVFuuf1fv1G4airdPoflyTO7SIWEFeeiLlgy/xHWGeSDCkWLHOz9nHlVI7xsm/V/lybVI3Sin2NIe78xsAE6pd7OojIh9tCzCm3EFNRf/VVk1tiRPriUJZI4meTr69PcuOLgPDgIoyfQnLNfo/YAVqOHoiURHJpzxOBqzfFiAQUhwyu5DD5hTx/ib/kMOM7KLda+APqF6tQCbWuPp5Iht3hHol1WPUVJqeSN9GjC2jQEQqShM3YWyL5kgqS0f26x8OaBGxAiMfPZEyiAy2OmtkhbNiLdPnTC/gwGkevD7VXdGca/Z295Pq8UTGVbtobI10FxEGggY79pjNBPtSXeEkGFJ0dPX+Np6ob9ZIw+MWSoukXzirNdoKpaJUX8Jyjf4PWEG3J5JnIjKKEutbd4coLzFzDbHE9FCzOOwi1uIkXkTGVjoxjJ6E8ZZdIQwFM+v6i0hsGWt8XsQfNPD6Vb9q9ZFIZZmzX+uTWO1Ipc6J5Jy8fgeKyAwR+b2IPJVrWwal2xMZhuGsEbLEd+vuENMmuBERJo9z43DAtj2Je07ZTVNUKKrjch2xBoWNUWGI1ZHsl8ATqYn2worPi7QMUK0+Eqkq79/6pC3qiVTqnEjOsf0/ICIPikiDiKzpc/8pIrJeRDaJyA0ASqnNSqkr7LYxZfLVEyEExgAdYEfY6qwde8NMjS6ddbuEiTWu7ul+uSYWv48PvcS8i8boeNtNO4KUFAoTqvuLQkx84j2RwWpERhpVCVqftMbOacnIf/35Ti5k/CHglPg7RMQJ3AecCswBLhCROfablib56IkM1Q5+BHkiXp9Bu9fotfppyng32/fkR06ktcOs54ivHq/pE6LaVB9kRl3i4UpaRJz9Vme1dUQoKRzZNTLDBdtFRCn1KtDc5+7DgE1RzyMIPA6cZbdtaZO3nggDL/M1OkE85s8wZ2805xD/LX7SWBd7msKDjpa1izav0b3KKEZFqQO3CxqjTQTrG8JMG5+4pbnHLVSUOkatiFSWOWn3GoQjPf/L1k6DCp0PyQvyJaA4CdgR93c9MElEqkXkfuAQEblxoCeLyJUiskJEVjQ2Nmbb1v7koycyVCdfo3PELO/dHe1wO666xxOprTJ7LsVmfueSto5Iv1VEIkJNpYt9bRE6ukxPalLtwN1vqyuc7IvLiTS3RRAZHauTYrUibXGtT9o6DSpHwWsfDuT1f0Ep1aSU+qpSaqZS6o5BtntAKbVIKbVo7NixdpoYNSCfPZEBRESNnA6+DQmW0MbakTe0JG4jbietnUbCVURjK500tkTY2WDmbiaNHVhEzFqRnvDc3uYINRVOnMO8wWIyxM5dfOuT1s5IP+9OkxvyRUR2ApPj/q6L3jc8MHwgBSD5cjpJzhMZISLS3B7B4aDXN9OYoAw0QdBOBvrWPLbKSWNrpLueZeJQIhLniexuCo+aYUy10SFb8RX+rR3aE8kX8uW/8C4wS0Smi4gHOB94Psc2JY/qyi8vBJIUkZERzmppj1BZ6sDh6PlWHrvwxAr9coVSirYBvjWPjXoX3SIyiCjUVJhzNWJ5gd37Ro+IzJjkweXsGdjVfU51TiQvyMUS38eAN4H9RaReRK5QSoWBa4C/AeuAJ5VSa+22LW0MX37lQ6DHyzASj2EdSQOpmhNUbleUOvC4pXv6Xa7w+hXhSOLcRU2li1AYPtwSYGylkwLPwB/HmkoXSpm5kEDQYF9rZNSIiMctzKzzsG6LuVzd6zPPaZWuEckLbH8XKqUuGOD+ZcAym82xhnz0RGJJ89hS3r4YneCybhZ2LmnpMPqtUhIRaquc3fmSXNEWrW9IFHqJLfN9f2OAA6YNvkou1pSxoSXClt1mDiW+ZfxIZ850D39500vEUHEtT7Qnkg9oKbeCfPZE1ACeiBpZ4axELcFrq5w590S6i+ISJdajITd/UA2aVAeoG2c+vmNviJfe9VLoEQ6ZXWixtfnLgdMK8AcUW3eFelqe6JxIXqD/C1aQl55IAeAc3BMZAeEspRQtHZGEoY3aMa6c50Ri35oTJtbjBipNnzi4VzGhxoXHLTz5UgcvvdvFGceUjqpCuzkzCgD44ONAXMsT7YnkA1pErCAfPRER09MY4TkRr08RCsOYiv4XlJpKJy3tkZy2hG/vbnnS3774PM60iYkLDWM4HcLkcS627Q4xs87NZadVWGtonjOh2kltlZNVG/w93p32RPIC/V+wgnz0RMAUiUSeiFIjxhNpjuYcEoWzKkudGIp+LdTtZLALntMhHDXPfN/MnDS4iABcdloFnz22lDu+VktR4ej66IoI82cV8MGmwKB5Jo39jI7lHdnG8IE7zzwRMJPrCUUkAERGRE6kJVo7kWiuRuzC3drRv+2IXbR2RChwC0UFiS94t32lhpb25Ow7an4xR83Pw/eZTRwwrYB/vNPF6k0BKkodg65m09iH/i9YQT57IokS67H7ZPh7Ii3dnkj/t3IsZh7LS+SCtk5j0BGuIpIwFKfpzwHR1Wjvfuhn2oShPTeNPWgRsYJ8zIlANCeSwBMZQQOpmgeZq1EZ54mkw+594X4jbFOlrTNCRYn+mFnB9EluYk2Opw7QrFJjP/rdbQX57IkkSqyPoDbwbZ1mI8KyBBfqTDyRxtYwF92yi8tu20W7N31PZqC+WZrUKfQ4iDVlnpFEDkljD1pErEB7Ijmj02dQUig4Hf2Xu8ZyIvHdX5PlH2+Z4hsKwz/eHmCFWxKYLU/0x8xqFs3Jwy9toxT97s4UFQZC+emJyECeyMiZatjZZVBanPht7HIKpUXSb7RqMrzxgY8DpnmYPtHNGx/40ravtdOgUldWW8aNl1Zz3MLiQfuMaexF/ycyJTZLRPLUE1GDeSLDP5zV0WVQWjTwd6HKMmfKnkhnl8FHW4NcfGo5gaBiySsdBEMq5eK+QNDAH1B6KaqFnHR4CScdPvzftyMJ/e7OlNgsEUceeiJD5URGwOosr08N6ImAKSJ953MPxeZdQQwFB04vYPI4N+EIKe8DesJo5doT0YxgtIhkSt57In5QfS6AagSFs3yDeyIVpY6UPZHYbPap493dy3PTCYnFCg21J6IZyeh3d6bkuycC/b2REZRY7xgkJwLmBTzV1Vlbd4coLDC7AMfyGa1pJOd1jyfNaECLSKaoPPZEBmoHP4JyIp0+g7JBRcTMiaTSP2vb7hBTx7lxOKTbi2hLxxPp0D2eNCMf/e7OFGMYeCJ9q9aNzug43+G9riIYUvgDg+dECgsEw4BwChqwbU+IKePNc9M93zuNgsU2b2zuhf6YaUYu+t2dKfnsiTgG8kRGRgffHXvN4UyDzeJItVl6p8+cGjg12lajuFBwu3pCU6ng9Znez2A5G41muKPf3ZkyHDyRROEsGf6hrO17TBGxsgXG3qbovPOx5j5FhOJCB13+9NvJy+gZ+6EZhWgRyZRh4Yn0CWepkdEGftueECJQV2tdWC7WVsMV98nQGqDRDIwWkUzJZ09EBvJERkY4yxdQFHpEtwTXaHKI/vRlynDwRPol1vO0YaRGoxl2aBHJlLz2RArNW8Pf5wEFov/1Go0mc/SVJFO6PZF8FBEdzddoNNlFi0imGF3Rmgt9KjUazehDX/kyRfnyMx+i0Wg0NqBFJFOMrvzMh2g0Go0NaBHJmAh6LItGoxmtaBHRaDQaTdpoEdFoNBpN2mgR0Wg0Gk3aaBHRaDQaTdpoEdFoNBpN2mgR0Wg0Gk3aaBHRaDQaTdpoEdFoNBpN2mgR0Wg0Gk3aiFLpj/3MR0SkA1ifazv6UAPsy7URfdA2JU8+2qVtSg5tU/Lsr5QqS/VJI7Ffx3ql1KJcGxGPiKzQNg1NPtoE+WmXtik5tE3JIyIr0nmeDmdpNBqNJm20iGg0Go0mbUaiiDyQawMSoG1Kjny0CfLTLm1Tcmibkictu0ZcYl2j0Wg09jESPRGNRqPR2MSwFREROUVE1ovIJhG5IcHjBSLyRPTxt0VkWh7YdJmINIrIqujPl7Jsz4Mi0iAiawZ4XETk3qi974vIwmzak4Jdx4lIW9x5uiXL9kwWkVdE5EMRWSsi1yXYxvZzlaRddp+rQhF5R0RWR236vwTb2PrZS9ImWz97ccd1ish/RGRpgsdsv0YlYVPq50kpNex+ACfwMTAD8ACrgTl9tvkacH/09/OBJ/LApsuAX9p4nj4FLATWDPD4Z4C/AAJ8Eng7T+w6Dlhq43maACyM/l4GbEjwv7P9XCVpl93nSoDS6O9u4G3gk322sfuzl4xNtn724o77DeDRRP8ju89TkjalfJ6GqydyGLBJKbVZKRUEHgfO6rPNWcAfo78/BZwgIpJjm2xFKfUq0DzIJmcBf1ImbwGVIjIhD+yyFaXUbqXUyujvHcA6YFKfzWw/V0naZSvR198Z/dMd/embWLX1s5ekTbYjInXAacDvBtjE7mtUMjalzHAVkUnAjri/6+n/4ereRikVBtqA6hzbBHBuNBzylIhMzqI9yZCszbngiGh44i8icpBdB42GFA7B/DYbT07P1SB2gc3nKhoOWQU0AP9QSg14rmz67CVjE9j/2bsH+F/AGOBx289TEjZBiudpuIrIcOUFYJpSah7wD3q+hWh6sxKYqpSaD/wCeNaOg4pIKbAEuF4p1W7HMZNhCLtsP1dKqYhSagFQBxwmInOzfcyhSMImWz97InI60KCUei+bx0mFJG1K+TwNVxHZCcQrZF30voTbiIgLqACacmmTUqpJKRWI/vk74NAs2pMMyZxH21FKtcfCE0qpZYBbRGqyeUwRcWNeqB9RSj2dYJOcnKuh7MrFuYo7divwCnBKn4fs/uwNaVMOPntHAWeKyFbM0PbxIvLnPtvYfZ6GtCmd8zRcReRdYJaITBcRD2ZS6vk+2zwPXBr9/b+Al1U0c5Qrm/rE0M/EjHHnkueBS6Irjz4JtCmldufYJkRkfCw2LCKHYb5Ps/bhih7r98A6pdTPBtjM9nOVjF05OFdjRaQy+nsRcBLwUZ/NbP3sJWOT3Z89pdSNSqk6pdQ0zGvBy0qpi/tsZut5SsamdM7TsGzAqJQKi8g1wN8wV0U9qJRaKyK3ASuUUs9jfvgeFpFNmEnc8/PApmtF5EwgHLXpsmzaJCKPYa7eqRGReuBWzKQjSqn7gWWYq442AV3A5dm0JwW7/gu4SkTCgA84P8tfAI4CvgB8EI2rA3wHmBJnUy7OVTJ22X2uJgB/FBEnpmA9qZRamsvPXpI22frZG4gcn6dkbEr5POmKdY1Go9GkzXANZ2k0Go0mD9AiotFoNJq00SKi0Wg0mrTRIqLRaDSatNEiotFoNJq00SKi0Wg0mrTRIqLRaDSatNEiotEkiYi8KiLtIvKRiJyR4nOLRORf0YK4gbaZKSIf9LmvQES2iMj86PGHZYGwZuSiRUSjSZ6ngLuUUgcopV5I8blfBJ5WSkUG2WYLUCci8Z/LK4FXlVKrgZeAz6d4XI0mq2gR0WiSZy4wUUReFpEdInIBgIicKSJL4jcUkatE5Bdxd10EPBf3+HQReU5EVog5lW9/pZQBbAemRbcpAr6J2RYGzA69F2XptWk0aaFFRKNJnrlAh1LqeEzP4gvR+2+n50If42PgQIBoQ84ZSqmt0b/dmB1Sv6GUWgR8D4iNU14HHBD9/WrghdjzgDXAJyx9RRpNhuj4qkaTPLOAk6O/+wCviMwHHEqpNSIyFfiMUurX9J6uVwO0xu3ns8BBwJJoA14X8Fr0sXXA/iLyKnANcHjsSUqpiIgERaQsOulQo8k5WkQ0miQQkSnAzriL91xMz2ABEBvycxKm0ADMAVZHf/cBhXG7mw98Vyn1+wSHWgecAFyHOUNkb5/HCwB/+q9Eo7EWHc7SaJJjLj2iAKYQvI/5GSqNrro6ByiL5jIuAx4FUEq1AE4RiQnJbuDkWAJdRA6OzQTBFJHDMMNlP443QESqgX1KqZD1L0+jSQ8tIhpNcszFFI0Y86J/LwNmAKuA+zHDVCuAB5RSK+O2/ztwdPT3BzE/e+uic0K+HTcDZANwcPT5rX1sWAy8aM3L0WisQc8T0WhsQEQWAl9XSn1hyI0H3sfTwA1KqQ3WWabRZIb2RDQaG4h6Ja8MVmw4GNEVXs9qAdHkG9oT0Wg0Gk3aaE9Eo9FoNGmjRUSj0Wg0aaNFRKPRaDRpo0VEo9FoNGmjRUSj0Wg0aaNFRKPRaDRpo0VEo9FoNGnz/wHPu7msdqFhLwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Plot the absorption coefficient for the both systems\n", "import matplotlib.pyplot as plt\n", "\n", "silicon_equil.abs_calc\n", "data = np.loadtxt(str(os.getcwd())+'/'+f\"si_zg/eps/abs.dat\")\n", "data2 = np.loadtxt(str(os.getcwd())+'/'+f\"si_equil/eps/abs.dat\")\n", "plt.plot(data2[:, 0], data2[:, 1],color='royalblue', label='Equilibrium structure')\n", "plt.plot(data[:, 0], data[:, 1],color='gold', label='ZG structure')\n", "plt.xlim(0.0,4.5)\n", "#plt.ylim(1e-3,1e2)\n", "plt.yscale('log')\n", "plt.legend()\n", "plt.xlabel(r'$\\hbar\\omega (eV)$')\n", "plt.ylabel(r'$\\alpha(cm^{-1})$')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\underline{\\large{\\rm Question}}$\n", "\n", "Can you observe differences in the optical absorption? As an exercise, (i) investigate convergence of the spectra with respect to the k-point grid and (ii) investigate the temperature dependence of the spectra." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.7" } }, "nbformat": 4, "nbformat_minor": 4 }