estudio ab-initio de las propiedades estructurales,

elásticas y piezoeléctricas del compuesto CdSe

ESTuDIO AB-INITIO
ab-initio DE LAS
study PROPIEDADES elastic
of structural, ESTRuCTuRALES,
and
ELáSTICAS y PIEzOELéCTRICAS DEL COMPuESTO
piezoelectric properties of CdSe compound CdSe
† †‡
N. V. González Lemus y E. A. Albanesi
†
INTEC-CONICET, GÜEMES 3450, 3000 SANTA FE, ARGENTINA,
‡
FACULTAD DE INGENIERÍA, UNIVERSIDAD NACIONAL DE ENTRE RÍOS, 3101 ORO VERDE, ARGENTINA
e-mail: vanessagonzalez@santafe-conicet.gov.ar

Recibido 30/03/2012; aprobado 10/09/2012

En el presente trabajo presentamos un estudio teórico de la estructura electrónica, propiedades ópticas, elásticas y
piezoeléctricas del semiconductor con estructura wurtzita CdSe, obtenidas mediante cálculos ab-initio empleando
el método de ondas planas y pseudopotenciales dentro de la teoría de la funcional densidad (DFT), utilizando la
aproximación de gradiente generalizado(GGA). los cálculos basados en DFT suelen despreciar el valor del gap de
energía en muchos sistemas, para corregir este efecto se introdujo la corrección GW en el marco de la teoría de
perturbación de muchos cuerpos, a partir de esto calculamos el modulo de bulk ajustando valores de energía y
volumen mediante la implementación de la ecuación de Birch-Murnaghan, por ultimo mediante las teorías de
Voigt y Reuss, estimamos el módulo de elasticidad transversal, el módulo de Young, el coeficiente de Poisson y
la corrección anisotrópica del módulo de Bulk.

Palabras Claves: CdSe, ab-initio, DFT, GW, Elasticidad.

In this paper we present a theoretical study of electronic structure, optical, elastic and piezoelectric properties of
semiconductor with wurtzite structure CdSe, obtained by ab-initio calculations using the method of plane waves
and pseudopotentials within the density functional theory (DFT ) using the generalized gradient approximation
(GGA). DFT calculations frequently disparage the value of the energy gap in many systems, to correct this, GW
correction was introduced in the framework of perturbation theory of many bodies, from this we calculate the
bulk modulus adjusted the energy and volume values through the implementation of the Birch-Murnaghan
equation, finally by the theories of Voigt and Reuss, we estimate the elastic modulus, Young's modulus, Poisson's
ratio and anisotropic correction of bulk module.

Key Word: CdSe, ab-initio, DFT, GW, Elasticity

I. INTRODuCCIÓN las propiedades ópticas la función dieléctrica (parte real

El seleniuro de cadmio (CdSe) es un compuesto y compleja), índice de refracción , reflectividad,
semiconductor del grupo II-VI, con estructura wurtzita, coeficiente de absorción y conductividad; e
grupo espacial p63mc Fig 1. Este compuesto ha cobrado implementamos el cálculo perturbativo GW(4), aplicado
gran interés en los últimos años debido a sus múltiples e a la corrección de la estructura de bandas y obteniendo
importantes aplicaciones tales como los detectores de una corrección del valor del band gap.
radiación nuclear que operan a temperatura ambiente,
dispositivos ópticos no lineales, celdas solares,
transistores, diodos emisores de luz, dispositivos
electroluminiscentes, etc.
En la actualidad conocer la estructura electrónica de
los materiales mediante cálculos a primeros principios o
ab-initio, permite la comprensión de la naturaleza de
los compuestos así como su comportamiento a escala
cuántica y predecir su respuesta a escalas mayores,
debido a que muchos de estos parámetros están
relacionados con fenómenos ópticos, magnéticos,
eléctricos, estructurales, etc.; por lo tanto, la teoría de
elasticidad del medio continuo, puede ser abordada
también cuánticamente, ya que desde el punto de vista
microscópico, la elasticidad tiene su origen en las
fuerzas intermoleculares.
Presentamos en este trabajo el cálculo ab-initio de la
estructura electrónica para el CdSe, basado en la teoría
de la funcional densidad (DFT) (1) con la aproximación
GGA (2) para el cálculo del término de correlación e
intercambio, a través del método de ondas planas y Figura 1. Vista superior del CdSe estructura wurtzita,
pseudopotenciales(3), también presentamos un estudio de parámetros de red a=4,302Åy c=7,021 Å.

ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011 65

 Finalmente mostramos el estudio de las propiedades *
2

1  N (w )
elásticas como funciones de respuesta (RF) 
R (w ) 
perturbativas y mediante las teorías de Voigt y Reuss y *
1  N (w ) (5)
cálculo del modulo de bulk (5), y su corrección 
anisotrópica.
II. MéTODO DE CáLCuLO 2w 1 ( w)  û ( w) û1/2 2w
 ( w)  ( )  k ( w) (6)
Estructura electrónica y corrección GW c 2 c
Por otra parte una de las propiedades mecánicas más
En el marco de la DFT(1), utilizamos la utilizadas para el cálculo asociado a materiales es el
aproximación de gradiente generalizado GGA (2) para el módulo de bulk, que muestra principalmente la
termino de intercambio y correlación y la aproximación capacidad de restitución de un sólido sometido a presión
del pseudopotencial (3), incluida en el software hidrostática, también conocido como modulo de
ABINIT(6), donde se toma en cuenta que los electrones incompresibilidad, por definición el modulo de bulk está
de valencia determinan mayormente las propiedades dado por:
químicas y físicas de los sistemas; debido a que el V
B  V
esquema general de cálculo de la DFT puede subvalorar P (7)
el valor del gap de energía en muchos sistema, para Donde P es la presión y V el volumen del sistema, el
corregir este efecto se utilizo la corrección GW inverso del modulo de bulk corresponde a la
desarrollada por Hedin en el marco de la teoría de compresibilidad del material. Una ecuación de estado
perturbación de muchos cuerpos en la que se aproxima fue derivada por Murnaghan , y mejorada junto a birch
el valor de autoenergía utilizando un término que en 1947(10,11).Esta escribe la energía en función del
incluye la función de Green que describe la propagación volumen del sólido, y permite ajustar los resultados
de una partícula en el sistema de muchos electrones y el obtenidos mediante el cálculo de la energía del sistema a
potencial coulombiano apantallado(7). Para el cálculo diferentes variaciones de volumen para así obtener el
auto-consistente necesario para obtener la estructura valor del modulo de bulk(12), esta ecuación está dada
electrónica se utilizaron 112 puntos k en la 1BZ por:
generados de una grilla de 10X10X14. 9V 0 B 0
E (V )  E 0 
16
Propiedades ópticas y modulo de bulk
  V 2 / 3 3
  V 2 /3  
2
V 2 /3  
Las funciones dieléctricas se calcularon con las    0  1  B 0'   0  1   6  4 0   (8)
V  V V
autofunciones correspondientes al cálculo       

autoconsistente, calculando el elemento de matriz Donde B0' es un parámetro de ajuste y B0 es el

correspondiente al operador momento dipolar que para modulo de bulk. Para esto calculamos la energía para
la parte imaginaria del tensor dieléctrico complejo (8) diferentes variaciones de volumen desde -6% hasta
está dado por: +6%.
4π 2 e 2 Propiedades elásticas y piezoeléctricas
 2 ( w )   2 2    f ∣ p ∣ i  f ∣ p  ∣ i 
m w i, f Para el cálculo de las matrices elásticas, se relajo la
(1) estructura, haciendo que la fuerza sobre los núcleos sea
 W i (1  W f ) ( E f  E i    ) d 3 k
nula, configurando la celda en equilibrio, esta fuerza se
Donde <f / pj / i> es la matriz dipolar, y f, i son los
halla mediante el teorema de Hellmann-Feymann (13,14),
estados final e inicial respectivamente. Wi es la función
dado por:
de distribución de Fermi para el estado i- ésimo, y Ei es dEe
la energía del electrón en el i-ésimo estado, siendo α y β Fi  
dRi (9)
componentes en coordenadas cartesianas, entonces el
tensor dieléctrico complejo es de la forma: Y el estrés de la estructura se obtiene mediante el
teorema de Nielsen-Martin (15):
 (w)  1  i 2 ( w)  N 2 (2)
1  p  p i ,β 
Donde las partes real e imaginaria son funciones del σ αβ    i,α   ri,α  i,α Vij 
índice de refracción complejo N. De la parte imaginaria Ω i  mi j i  (10)
del tensor dieléctrico se obtiene la correspondiente parte Donde E es la energía del sistema R son las
real por medio de las relaciones de dispersión de coordenadas nucleares, P es el momento lineal y V es el
Kramers-Kronig: termino de interacción colombiano. La función de onda
2  w '  2 ( w ') de primer orden se calcula para cada átomo en todas las
1 ( w)  1  P  dw ' (3) direcciones; esta función de onda se utiliza para calcular
 0 w '2  w2
El índice de refracción complejo esta dado por: la segunda derivada de la energía E, con respecto a
N ( w)  n( w)  ik ( w) (4) todos los pares de desplazamientos atómicos fuera de
sus posiciones de equilibrio, este es el tensor de fuerzas
Donde la parte real es el índice de refracción y la
interatómicas K:
parte imaginaria representa el coeficiente de extinción
(9)
. La reflectividad R y el coeficiente de absorción α d 2E
K ij   0
están dados por: dui du j
 ,
(11)

66 ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011

 Donde Ω0 es el volumen de la celda, Ɛ el campo
eléctrico homogéneo y η la deformación homogénea.
Posteriormente se calculan las segundas derivadas de E
con respecto al desplazamiento de un átomo y con
respecto a uno de los componentes de ε este forma el
tensor dinámico carga efectiva (carga efectiva de Born):
d 2E
Z ia    0
du i d  a  (12)
Donde a, representa la coordenada cartesiana de la
componente. Este tensor se utiliza para calcular la
contribución atómica a la relajación los tensores
elásticos y piezoeléctricos (16), se computan las
funciones de onda para una perturbación de
deformación aplicada sobre el sistema que se utilizan
para calcular tres tipos de segundas derivadas de energía
E, primero las derivadas con respecto a dos
componentes de la tensión, que dan como resultado el
tensor elástico de iones congelados: Figura 2.Estructura de bandas del CdSe calculada con
la aproximación GGA.
d 2E
C 
d   d   u , (13) De la DOS se puede notar que el grupo de bandas
más profundas corresponden a los orbitales s de Se
Donde α y β=1,2…6 en notación de Voigt. En principalmente, con una leve contribución de los d del
segundo lugar se calculan las derivadas de E respecto a Cd Fig.3, el grupo de bandas que inicia en -8,5 eV está
un esfuerzo de deformación y una componente del compuesta casi totalmente por los orbitales d del Cd el
campo eléctrico, para hallar el tensor piezoeléctrico de grupo superior de bandas está conformado por en su
iones congelados (frozen-ion piezoelectric tensor): mayoría por los orbitales p del Se.
d 2E
ea  (14)
d  d
a  u
Por último se obtienen las segundas derivadas de E
con respecto a un esfuerzo y un desplazamiento
atómico, que da como resultado el tensor de fuerzas
internas (17):
d 2E
 i    0
du i d   (15)


III. DISCuSION DE RESuLTADOS

Estructura electrónica y corrección GW
La estructura electrónica del CdSe se calculo
empleando pseudopotenciales Troullier-Martins, con
una energía de corte de 20 eV para el cálculo
autoconsistente GS. Los cálculos de la DOS abarcan en
un rango de energías de -10 ev a 10 eV a lo largo del
camino Γ - A - H – K - Γ - M - L - A en la zona de Figura 3.Densidad de estados del CdSe. Comparación
Brillouin (1BZ). de la densidad de estados del Cd, Se y DOS total del
En la estructura de bandas Fig. 2 se puede observar compuesto.
que el sistema es un semiconductor de gap directo, con
un Eg de 0,41eV, el gap se encuentra en el punto La estructura electrónica reportada a partir de
Γ(0,0,0). La estructura de bandas presenta cálculos GGA con la aproximación de
degeneraciones en las bandas en los caminos en la 1BZ: pseudopotenciales subestima el valor del band gap
Γ - A y H - K - Γ – L - M, donde las bandas están Eg=0,408 eV, usando la aproximación GW de la teoría
doblemente degeneradas, se pueden diferenciar 2 grupos de perturbaciones, se calculo la corrección para las
principales de bandas, la primera en el rango de energía primeras 4 bandas de conducción tomando un rango de
que va desde -8,5 hasta -7,3 eV y la segunda en el rango energía de -4,5 a 7,5 eV en el camino K - Γ- M Fig. 4,
de -4,8 a 0 eV. obteniéndose un valor de gap Eg=1,67eV más próximo
a los valores experimentales reportados.

ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011 67

 En la parte real se observan 2 picos, el pico más
pronunciado se encuentra a 3,8 eV y el menor a 5,8 eV
en las direcciones xx e yy, mientras que en la dirección
zz se encuentras 2 picos levemente corridos, uno a 4 eV
y otro a 5,5 eV. En la parte imaginaria se encuentran 2
picos notables en todas las direcciones el primero en 4,2
en la dirección zz y en 4,4 en las direcciones xx e yy, y
el segundo en 6,3 para xx, yy y zz.
Tabla 2. Función dieléctrica presentada en este trabajo en
comparación con otros autores.
pico A1 A3 A2 A 4 B 1 B3 B2 B4
presente trabajo 4,3 4,3 6,2 6,3 4,7 4,5 6,4 6,4
Experimental(a) 4,2 5
Excitonic(a) 4,4 5,7
Figura 4. Comparación de la estructura de bandas en el Experimental(a) 4,3 8,7
camino K-Γ-M obtenida con DFT y con la corrección GGA(b) 4 6,3 5,5 10
GW. LDA(c) 5,5 5,3
(a)
Ref(20)(b)Ref(21)(c)Ref(22)
Tabla 1.Band gap calculado en comparación con los
obtenidos experimentalmente por otros autores. Los picos representados en la grafica anterior de la
(a) (b) parte real e imaginaria de la función dieléctrica que
Presente trabajo Experimental Experimental
incluye la corrección obtenida mediante la GW
Band gap(eV) 1,67 1,74 1,67
(a) presentan un buen ajuste respecto a los obtenidos
Ref (18)(b)Ref(19)
experimental y teóricamente por otros autores y tienen
Propiedades ópticas y modulo de bulk correspondencia con la grafica de índice de refracción y
coeficiente de extinción Fig. 5 respectivamente, los
Las propiedades ópticas del CdSe se hallaron
picos de la curva de la parte imaginaria se ven
analizando la convergencia empleando 84 y 112 puntos
correspondidos además en la refelctividad Fig. 6 debido
k, hallando que se obtiene convergencia con 84 puntos k
a su relación matemática; se puede observar también
en la 1ZB sin la inclusión de la interacción spin-orbita.
que el material se comienza a comportar casi
Se calculo entonces la función dieléctrica, donde para se
isotrópicamente después de los 20 eV.
puede notar que el CdSe presenta isotropía en dirección
La representación de la conductividad óptica y el
xx e yy en todo el rango de energía estudiado tanto en la
coeficiente de absorción guardan una estrecha
parte real como imaginaria de la función Fig. 5.
coincidencia en el comportamiento y distribución de
picos en el rango energético, presentando en todo en los
dos casos picos paralelos con diferentes intensidades en
las 3 direcciones.

Figura 6. (Niveles de arriba abajo)Coeficiente de

Figura 5.propiedades ópticas del CdSe. (Niveles de arriba absorción. Conductividad óptica Reflectividad.
abajo)Parte real e imaginaria de la función dieléctrica.
Índice de refracción. Coeficiente de extinción.

68 ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011

 Se halló un ajuste entre el los resultados de varias Donde Sc es el compliance tensor en unidades 10-
2
simulaciones, mediante métodos ab-intio basados en la GPa-1. En este caso las matrices elásticas Cc y Sc
DFT, en las que hallamos diferentes valores de volumen elásticos muestran un comportamiento anisotrópico del
y energía (calculamos la energía para diferentes tipo transversalmente isotrópico, debido a que las
variaciones de volumen desde -6% hasta +6%), estos constantes del tensor elástico cumplen con la relación
datos se fitearon por medio de la ecuación de Birch- C11=C22, C12=C21, C13=C23=C31=C32, C44=C55,
Murnaghan Fig. 7, obteniéndose un modulo de bulk de C66=(C11=C12)/2 y también que C33, es diferente y las
47 GPa. demás constantes tienden a 0.
A partir de la condición de isotropía transversal
pueden hallarse constantes elásticas: modulo de Young
en dirección Ex=Ey≠Ez, la relación de Poisson
νxz=vyz≠vxy y el modulo de cizalladura Gxz=Gyz≠Gxy, del
compliance tensor.
Tabla 3.cosntantes elásticas del CdSe.
E1=E2 E3 V13=V23 V12 G13=G23 G12
414,659 111,763 0,144 0,319 22,851 157,075

La matriz piezoeléctrica anisotrópica en unidades

Figura 7.Ajuste energía vs volumen mediante la C/m2 obtenida es:
ecuación de Birch-Murnagham.  0, 000 0, 000 0, 35800518 
  0, 000 0, 000 0, 35800935 
Propiedades elásticas y piezoeléctricas 
 0, 000 0, 000  0, 74612283 
ec   
Se realizo el cálculo perturbativo RF para los  0, 000 0, 38550299 0, 00000034 
 0, 38550827 0, 000 0, 00000109 
tensores piezo-elásticos del CdSe optimizando de la  
estructura mediante el algoritmo BFGS y mediante la  0, 000 0, 000  0, 00000265 

herramienta ANADDB (Analysis of Derivative Calculamos el modulo de bulk anisotrópico a lo

DataBase) obtuvimos el tensor elástico inverso largo de los ejes (a, b y c) a partir de sus constantes
(compilance tensor) que se obtiene por inversión del elásticas de acuerdo con las ecuaciones:
tensor elástico C   , el tensor de tensión interna que dP  dP Ba dP Ba
Ba  a  , Bb  b  , Bc  c 
representa el desplazamiento en respuesta a una tensión da 1     db  dc 
y por último los tensores elásticos y piezoeléctricos con 
Brel 
iones anclados y relajados. (1     )2 (16)
El tensor elástico con átomos relajados resulta de la Donde
combinación de la segunda derivada para tensión   C 11  2 C 12  C 22 2  2 C 13   C 13  2  C 33  2  2 C 33 
interna y las constantes de la fuerza interatómica;  
 C 22  C 12  C 11  C 13   ( C 11  C 12 )( C 23  C 12 )
mientras que el tensor piezoeléctrico con átomos  C 22  C 12  C 33  C 13   ( C 12  C 23 )( C 13  C 23 ) (17)
relajados toma en cuenta las segundas derivadas para la
tensión interna, constantes de fuerza interatómica, y Con a =1, hallándose Ba=Bb=36,277 GPa,
cargas efectivas de Born. Las matrices en notación de Bc=144,365 GPa y un modulo relativo de Br=16,113
voigth que describen entonces el comportamiento GPa, de acuerdo con su estructura cristalina el CdSe
elástico del material son: presenta mayor dureza en el eje cristalino c mientras es
igual en a y b.
 4, 6705891 1,5290886 0, 24104690, 0000008 0, 0000018 0, 000001
 1,5290886
 4, 6705811 0, 24104740, 0000002 0, 0000000 0, 0000017  IV. CONCLuSIONES
 0, 2410484
Cc  
0, 2410504 1,1363779 0, 0000024 0, 0000000 0, 0000039 

Se calculo la estructura electrónica del CdSe
 0, 000000 0, 000000 0, 0000000, 2476656 0, 0000004 0, 00000  obteniendo se que el band gap directo hallado mediante
 0, 00000 0, 000000 0, 0000000 0, 0000004 0, 2476690 0, 00000 
  la aproximación GGA esta subestimado
 0,0000016 0, 000002 0, 0000004 0, 0000000 0, 0000000 1,5707520 
comparativamente con los resultados experimentales
reportados por otros autores, esto se debe que la forma
Donde Cc es el tensor elástico clamped ion en
del potencial de intercambio y correlación no describe
unidades 102GPa que comprueba que el CdSe cumple
correctamente los autovalores de los estados
con la condición para la estabilidad mecánica de Born
electrónicos, se logra una mejora significativa de estos
para cristales hexagonales, que exigen que C44 > 0, C11
resultados introduciendo la corrección GW, en la que se
>|C12|, C11 C33 >(C13)2y C33(C11 + C12) >2(C13)2
aproxima el valor de autoenergía mediante una función
 0, 2411616 0,077158 0,034788 0,0000004 0,000002 0, 0000003 
de Green.
 0,077158 0, 2411620 0,034788 0,000000 0,0000006 0,000000 De la DOS mostramos que la ultima las bandas de

 0,034788 0,034789 0,8947477 0,0000086 0,0000001 0,000002 valencia están conformadas principalmente por los
Sc   
 0,000000 0,0000000 0,0000000 4, 0377030 0,0000065 0,0000000  orbitales p del Se hibridizados con los d del Cd.
 0,0000000 0,0000000 0, 00000 0, 0000065 4,0376472 0,0000000 
 
 0, 0000002 0,0000002 0,0000001 0,0000000 0,000000 0,6366377 

ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011 69

 Las propiedades ópticas presentan una anisotropía
biaxial para energías menores a 20eV después de la cual
el material tiene un comportamiento isotrópico.
Por último se obtuvieron las propiedades piezo-
elásticas donde se observo el comportamiento
transversalmente anisotrópico encontrándose que el
material es fácilmente polarizable en la dirección z; el
cálculo de corrección anisotrópica del modulo de bulk
agrega una corrección del 22,8% en los ejes cristalinos a
y b y del 64,7% en el eje c.
Referencias
1- Hohenberg, P. and Kohn, W. Physical Review B, 136 864
(1964).
2- Perdew, P. Burke, K. and Wang, Y. Physical Review B,
54 16553 (1996).
3- Phillips, J.C. and Pandey, K.C. Physical Review Letters,
30 787 (1973).
4- Rehr, J.J. Physica Scripta, T115 19 (2005).
5- Cohen, M.L. Physical Review B, 32 12 (1985).
6- Gonze, X. Beuken, J.M. Caracas, R. Detraux, F. Fuchs, M.
Rignanese, G.M. Sindic, L. Verstraete, M. Zerah, G. Jollet, F.
Torrent, M. Roy, A. Mikami, M. Ghosez, P. Raty, J.Y. and
Allan, D.C. Computational Materials Science, 25 478 (2002).
7- Hedin, L. Physical Review, 139 796 (1965).
8- Cabuk, S. Akkus, H. Mamedov, A.M. Physica B, 394 81
(2007).
9- Suleyman, C. Sevket, S. Central European Journal of
Physics, 6(3) 730 (2008).
10- Murnaghan, F.D. Proceedings of the National Academy of
Sciences of the United States of America, 30 9 244 (1944).
11- Birch, F. Physical Review, 71 809 (1947).
12- Ziambaras, E. Schroder. E. Physical Review B, 68 064112
(2003).
13- Hellmann, H. Einführung in die Quantenchemie. Leipzig:
Franz Deuticke, 285 (1937).
14- Feynman, R. P. Physical Review, 56 (1939).
15- Nielsen, O.H and Martin, R. Physical Review Letters, 50
697 (1983).
16- Mousavi, M. Kompany, A. Roknabad, MR. Arabshahi, H.
Indian Journal of Science and Technology, 3 6 (2010).
17- Vanderbilt D. Notes relevant to the ANADDB module of
ABINIT, (2004)
18- Kucur, E. Riegler, J. Urban, G. and Nann, T. Journal of
Chemical Physics, 119 4 (2003).
19- Patidar, D. Rathore, S. Saxena, S. Sharma, T.
Chalcogenide Letters, 5 ( 2008).
20- Kim, D. Klein, S. Ren, F and Chang, H. Physical Review,
8(1994).
21- Huang, P, Cheng, D. Wu, S. Hu, M. Shen, J. Appl. Phys.
Lett. 90(2007).
22- Wang, L. and Zunger, A. Physical Review B, 53(1996).

70 ANALES AFA Vol. 23 N.2 (65-70) MONTEVIDEO 2011