El método de Monte Carlo para el transporte de electrones es un enfoque de Monte Carlo (MC) semiclásico para modelar el transporte de semiconductores . Suponiendo que el movimiento del portador consiste en vuelos libres interrumpidos por mecanismos de dispersión, se utiliza una computadora para simular las trayectorias de las partículas a medida que se mueven a través del dispositivo bajo la influencia de un campo eléctrico utilizando la mecánica clásica . Los eventos de dispersión y la duración del vuelo de las partículas se determinan mediante el uso de números aleatorios.
El modelo de ecuación de transporte de Boltzmann ha sido la principal herramienta utilizada en el análisis del transporte en semiconductores. La ecuación de BTE viene dada por [ cita requerida ] :
La función de distribución , f , es una función adimensional que se utiliza para extraer todos los observables de interés y da una descripción completa de la distribución de electrones tanto en el espacio real como en el espacio k . Además, representa físicamente la probabilidad de ocupación de partículas de energía k en la posición r y el tiempo t . Además, debido a que es una ecuación integro-diferencial de siete dimensiones (seis dimensiones en el espacio de fase y una en el tiempo), la solución a la BTE es engorrosa y se puede resolver en forma analítica cerrada bajo restricciones muy especiales. Numéricamente, la solución a la BTE se emplea utilizando un método determinista o un método estocástico. La solución del método determinista se basa en un método numérico basado en cuadrícula, como el enfoque de armónicos esféricos, mientras que el Monte Carlo es el enfoque estocástico utilizado para resolver la BTE.
El método semiclásico de Monte Carlo es un método estadístico utilizado para obtener una solución exacta a la ecuación de transporte de Boltzmann, que incluye una estructura de bandas compleja y procesos de dispersión . Este enfoque es semiclásico porque los mecanismos de dispersión se tratan de manera mecánica cuántica utilizando la regla de oro de Fermi , mientras que el transporte entre eventos de dispersión se trata utilizando la noción clásica de partícula. El modelo de Monte Carlo, en esencia, sigue la trayectoria de la partícula en cada vuelo libre y elige un mecanismo de dispersión correspondiente de manera estocástica. Dos de las grandes ventajas del método semiclásico de Monte Carlo son su capacidad para proporcionar un tratamiento mecánico cuántico preciso de varios mecanismos de dispersión distintos dentro de los términos de dispersión y la ausencia de suposiciones sobre la forma de distribución de portadores en energía o espacio k. La ecuación semiclásica que describe el movimiento de un electrón es
donde F es el campo eléctrico, E(k) es la relación de dispersión de energía y k es el vector de onda de momento. Para resolver la ecuación anterior, se necesita un conocimiento sólido de la estructura de banda (E(k)). La relación E(k) describe cómo se mueve la partícula dentro del dispositivo, además de representar información útil necesaria para el transporte, como la densidad de estados (DOS) y la velocidad de la partícula. Se puede obtener una relación E(K) de banda completa utilizando el método pseudopotencial semiempírico. [1]
Tanto los modelos de difusión por deriva (DD) como los hidrodinámicos (HD) pueden derivarse de los momentos de la ecuación de transporte de Boltzmann (BTE) utilizando una aproximación simplificada válida para dispositivos de canal largo. El esquema DD es el enfoque más clásico y generalmente resuelve la ecuación de Poisson y las ecuaciones de continuidad para portadores considerando los componentes de deriva y difusión. En este enfoque, se supone que el tiempo de tránsito de carga es muy grande en comparación con el tiempo de relajación de energía. [2] Por otro lado, el método HD resuelve el esquema DD con las ecuaciones de balance de energía obtenidas a partir de los momentos de BTE. [3] [4] Por lo tanto, se pueden capturar y calcular detalles físicos como el calentamiento del portador y el efecto de sobreimpulso de velocidad . No hace falta decir que se requiere un método de discretización preciso en la simulación HD, ya que las ecuaciones gobernantes están fuertemente acopladas y uno tiene que lidiar con un mayor número de variables en comparación con el esquema DD.
La precisión de los modelos semiclásicos se compara con base en el BTE investigando cómo tratan el problema clásico de sobreimpulso de velocidad, un efecto clave de canal corto (SCE) en estructuras de transistores. Esencialmente, el sobreimpulso de velocidad es un efecto no local de dispositivos escalados, que está relacionado con el aumento observado experimentalmente en el impulso de corriente y la transconductancia. [5] A medida que la longitud del canal se hace más pequeña, la velocidad ya no está saturada en la región de campo alto, sino que sobrepasa la velocidad de saturación predicha. La causa de este fenómeno es que el tiempo de tránsito de la portadora se vuelve comparable al tiempo de relajación de energía y, por lo tanto, las portadoras móviles no tienen tiempo suficiente para alcanzar el equilibrio con el campo eléctrico aplicado mediante la dispersión en los dispositivos de canal corto. [6] El resumen de los resultados de la simulación (Herramienta de Illinois: MOCA) con el modelo DD y HD se muestra en la figura de al lado. En la figura (a), se muestra el caso en el que el campo no es lo suficientemente alto como para causar el efecto de sobreimpulso de velocidad en toda la región del canal. Obsérvese que en dicho límite, los datos del modelo DD se ajustan bien al modelo MC en la región sin sobreimpulso, pero el modelo HD sobreestima la velocidad en esa región. El sobreimpulso de velocidad se observa solo cerca de la unión de drenaje en los datos MC y el modelo HD se ajusta bien en esa región. A partir de los datos MC, se puede notar que el efecto de sobreimpulso de velocidad es abrupto en la región de campo alto, que no está incluida correctamente en el modelo HD. Para condiciones de campo alto como se muestra en la figura (b), el efecto de sobreimpulso de velocidad se produce en casi todo el canal y los resultados HD y MC son muy similares en la región del canal.
La estructura de banda describe la relación entre la energía (E) y el vector de onda (k). La estructura de banda se utiliza para calcular el movimiento de los portadores bajo la acción del campo eléctrico, la tasa de dispersión y el estado final después de la colisión. La estructura de banda del silicio y su zona de Brillouin se muestran en la figura siguiente, pero no existe una expresión analítica que satisfaga toda la zona de Brillouin . Mediante el uso de alguna aproximación, existen dos modelos analíticos para la estructura de banda, a saber, el modo parabólico y el modo no parabólico.
Para el concepto de estructura de bandas, generalmente se suponen bandas de energía parabólicas por simplicidad. Los electrones residen, al menos cuando están cerca del equilibrio, cerca de los mínimos de la relación E(k). Entonces, la relación E(k) se puede extender en una serie de Taylor como
Debido a que la primera derivada se desvanece en el mínimo de la banda, el gradiente de E(k) es cero en k = 0. Por lo tanto,
lo que da la definición del tensor de masa efectivo
Esta expresión es válida para semiconductores que tienen una masa efectiva isotrópica, por ejemplo GaAs. En el caso del silicio, los mínimos de la banda de conducción no se encuentran en k = 0 y la masa efectiva depende de la orientación cristalográfica del mínimo como
donde describen la masa efectiva longitudinal y transversal, respectivamente.
Para campos aplicados más altos, los portadores se encuentran por encima del mínimo y la relación de dispersión, E(k), no satisface la expresión parabólica simple descrita anteriormente. Esta no parabolicidad se describe generalmente mediante
donde es un coeficiente de no parabolicidad dado por
donde es la masa del electrón en el vacío, y Eg es la brecha de energía. [7]
Para muchas aplicaciones, la estructura de banda no parabólica proporciona una aproximación razonable. Sin embargo, en caso de transporte de campo muy alto, se requiere un mejor modelo físico de la estructura de banda completa. Para el enfoque de banda completa, se utiliza una tabla de E(k) generada numéricamente. El enfoque de banda completa para la simulación de Monte Carlo fue utilizado por primera vez por Karl Hess en la Universidad de Illinois en Urbana-Champaign. Este enfoque se basa en el método pseudopotencial empírico sugerido por Cohen y Bergstresser [18]. El enfoque de banda completa es computacionalmente costoso, sin embargo, tras el avance de la potencia computacional, se puede utilizar como un enfoque más general. [8]
Para este tipo de simulación, se inyecta un portador y se sigue el movimiento en el dominio hasta que sale por contacto. Luego se inyecta otro portador y se repite el proceso para simular un conjunto de trayectorias. Este enfoque es útil principalmente para estudiar propiedades en masa, como la velocidad de deriva en estado estable como función del campo.
En lugar de una única portadora, se simula un gran conjunto de portadoras al mismo tiempo. Este procedimiento es obviamente un buen candidato para la supercomputación, ya que se puede aplicar paralelización y vectorización. Además, ahora es posible realizar promedios de conjunto directamente. Este enfoque es adecuado para simulaciones transitorias.
Este método combina el procedimiento de Monte Carlo de conjunto con la ecuación de Poisson y es el más adecuado para la simulación de dispositivos. Normalmente, la ecuación de Poisson se resuelve a intervalos fijos para actualizar el campo interno y reflejar la redistribución interna de carga debido al movimiento de los portadores.
La probabilidad de que el electrón sufra su próxima colisión durante dt alrededor de t está dada por
donde P[k(t)]dt es la probabilidad de que un electrón en el estado k sufra una colisión durante el tiempo dt. Debido a la complejidad de la integral en el exponente, no es práctico generar vuelos libres estocásticos con la distribución de la ecuación anterior. Para superar esta dificultad, la gente utiliza un esquema ficticio de “autodispersión”. Al hacer esto, la tasa de dispersión total, incluida esta autodispersión, es constante e igual a, digamos, . Por selección aleatoria, si se selecciona la autodispersión, k′ después de la colisión es igual a k y el portador continúa su vuelo sin perturbaciones. Al introducir una constante , la ecuación anterior se reduce a
Los números aleatorios r se pueden utilizar de forma muy sencilla para generar vuelos libres estocásticos, cuya duración estará dada por . El tiempo de computación utilizado para la autodispersión se compensa con creces con la simplificación del cálculo de la duración del vuelo libre. [9] Para mejorar la velocidad del cálculo del tiempo de vuelo libre, se utilizan varios esquemas como la “Técnica Constante” y la “Técnica por Partes” para minimizar los eventos de autodispersión.
Propiedades importantes de transporte de carga de los dispositivos semiconductores, como la desviación de la ley de Ohm y la saturación de la movilidad de los portadores, son una consecuencia directa de los mecanismos de dispersión. Por lo tanto, es de gran importancia para una simulación de dispositivos semiconductores capturar la física de tales mecanismos. La simulación Monte Carlo de semiconductores, en este ámbito, es una herramienta muy poderosa por la facilidad y la precisión con la que se puede incluir una gama casi exhaustiva de mecanismos de dispersión. La duración de los vuelos libres se determina a partir de las tasas de dispersión. Al final de cada vuelo, se debe elegir el mecanismo de dispersión apropiado para determinar la energía final del portador dispersado o, equivalentemente, su nuevo momento y ángulo de dispersión. En este sentido, se distinguirán dos grandes tipos de mecanismos de dispersión que se derivan naturalmente de la teoría cinética clásica de la colisión entre dos cuerpos:
Dispersión elástica , en la que la energía de la partícula se conserva después de ser dispersada. Por lo tanto, la dispersión elástica solo cambiará la dirección del momento de la partícula. La dispersión de impurezas y la dispersión superficial son, con una buena aproximación, dos buenos ejemplos de procesos de dispersión elástica.
Dispersión inelástica , en la que la energía se transfiere entre la partícula dispersada y el centro de dispersión. Las interacciones electrón-fonón son esencialmente inelásticas, ya que la partícula dispersa emite o absorbe un fonón de energía definida. Antes de caracterizar los mecanismos de dispersión con mayores detalles matemáticos, es importante señalar que, al ejecutar simulaciones de Monte Carlo de semiconductores, uno tiene que lidiar principalmente con los siguientes tipos de eventos de dispersión: [9]
Fonón acústico: el portador de carga intercambia energía con un modo acústico de vibración de los átomos en la red cristalina. Los fonones acústicos surgen principalmente de la excitación térmica de la red cristalina.
Óptica polar: El portador de carga intercambia energía con uno de los modos ópticos polares de la red cristalina. Estos modos no están presentes en los semiconductores covalentes. Los fonones ópticos surgen de la vibración entre sí de átomos de diferentes tipos cuando hay más de un átomo en la celda unitaria más pequeña, y generalmente son excitados por la luz.
Óptica no polar: la energía se intercambia con un modo óptico. Los fonones ópticos no polares generalmente deben considerarse en semiconductores covalentes y en el valle L de GaAs.
Fonón de intervalo equivalente: debido a la interacción con un fonón, el portador de carga pasa de estados iniciales a estados finales que pertenecen a valles diferentes pero equivalentes. Por lo general, este tipo de mecanismo de dispersión describe la transición de un electrón de un valle X a otro valle X, o de un valle L a otro valle L. [10]
Fonón Intervalle No Equivalente: Implica la transición de un portador de carga entre valles de diferentes tipos.
Fonón piezoeléctrico: Para bajas temperaturas.
Impureza ionizada: refleja la desviación de una partícula de su trayectoria balística debido a la interacción de Coulomb con una impureza ionizada en la red cristalina. Debido a que la masa de un electrón es relativamente pequeña en comparación con la de una impureza, la sección transversal de Coulomb disminuye rápidamente con la diferencia del módulo de momento entre el estado inicial y el final. [9] Por lo tanto, los eventos de dispersión de impurezas se consideran principalmente para la dispersión intravalle, la dispersión intrabanda y, en menor medida, la dispersión interbanda.
Portador-Portador: (interacciones electrón-electrón, hueco-hueco y electrón-hueco). Cuando la concentración de portadores es alta, este tipo de dispersión refleja la interacción electrostática entre portadores de carga. Este problema se vuelve muy rápidamente intensivo en términos computacionales con un número creciente de partículas en una simulación de conjunto. En este ámbito, los algoritmos Partícula-Partícula-Partícula-Malla (P3M), que distinguen la interacción de corto y largo alcance de una partícula con su gas de carga circundante, han demostrado ser eficientes al incluir la interacción portador-portador en la simulación de Monte Carlo de semiconductores. [11] Muy a menudo, la carga de los portadores se asigna a una cuadrícula utilizando un método de Nube en Celda, donde parte de la carga de una partícula dada se asigna a un número dado de puntos de cuadrícula más cercanos con un cierto factor de peso.
Plasmón: Refleja el efecto de la oscilación colectiva de los portadores de carga sobre una partícula dada.
Un enfoque computacionalmente eficiente para incluir la dispersión en la simulación de Monte Carlo consiste en almacenar las tasas de dispersión de los mecanismos individuales en tablas. Dadas las diferentes tasas de dispersión para un estado de partícula preciso, uno puede entonces seleccionar aleatoriamente el proceso de dispersión al final del vuelo libre. Estas tasas de dispersión se derivan muy a menudo utilizando la aproximación de Born , en la que un evento de dispersión es simplemente una transición entre dos estados de momento del portador involucrado. Como se discutió en la sección II-I, el problema cuántico de muchos cuerpos que surge de la interacción de un portador con su entorno circundante (fonones, electrones, huecos, plasmones, impurezas, ...) se puede reducir a un problema de dos cuerpos utilizando la aproximación de cuasipartícula, que separa el portador de interés del resto del cristal. [9] Dentro de estas aproximaciones, la Regla de Oro de Fermi da, al primer orden, la probabilidad de transición por unidad de tiempo para un mecanismo de dispersión de un estado a un estado :
donde H' es el hamiltoniano de perturbación que representa la colisión y E y E' son respectivamente las energías inicial y final del sistema constituido tanto por el portador como por el gas de electrones y fonones. La función de Dirac representa la conservación de la energía. Además, el término , generalmente denominado elemento de matriz, representa matemáticamente un producto interno de las funciones de onda inicial y final del portador: [12]
En una red cristalina, las funciones de onda y son simplemente ondas de Bloch . Cuando es posible, la expresión analítica de los elementos de la matriz se encuentra comúnmente mediante la expansión de Fourier del hamiltoniano H', como en el caso de la dispersión de impurezas [13] o la dispersión de fonones acústicos. [14] En el caso importante de una transición de un estado de energía E a un estado de energía E' debido a un fonón de vector de onda q y frecuencia , el cambio de energía y momento es:
donde R es un vector reticular recíproco . Los procesos Umklapp (o procesos U) cambian el momento de la partícula después de la dispersión y, por lo tanto, limitan la conducción en cristales semiconductores. Físicamente, los procesos U ocurren cuando el momento final de la partícula apunta fuera de la primera zona de Brillouin. Una vez que se conoce la probabilidad de dispersión por unidad de tiempo desde un estado k a un estado k', es interesante determinar la tasa de dispersión para un proceso de dispersión dado. La tasa de dispersión da la probabilidad por unidad de tiempo de dispersarse desde un estado k a cualquier otro estado en el espacio recíproco. Por lo tanto, la tasa de dispersión es
que se puede utilizar fácilmente para determinar el tiempo de vuelo libre y el proceso de dispersión, como se analiza en la sección 3-3. Es importante tener en cuenta que esta tasa de dispersión dependerá de la estructura de bandas del material (la dependencia surge de los elementos de la matriz).
Al final de un vuelo libre, se debe elegir aleatoriamente un modo y un ángulo de dispersión. Para determinar el mecanismo de dispersión, se deben considerar todas las tasas de dispersión de los mecanismos relevantes para la simulación, así como la tasa de dispersión total en el momento de la dispersión. La selección de un mecanismo de dispersión da como resultado simplemente la generación de un número aleatorio distribuido uniformemente 0 < r < 1 y la referencia a las siguientes reglas
Un enfoque computacionalmente eficiente para seleccionar el mecanismo de dispersión consiste en agregar un mecanismo de dispersión “vacío” de modo que permanezca constante en el tiempo. Si una partícula se dispersa de acuerdo con este mecanismo, mantendrá su trayectoria balística después de que se produzca la dispersión. Para elegir una nueva trayectoria, primero se debe derivar la energía (o momento ) de la partícula después de la dispersión.
donde el término representa la emisión o absorción de fonones y el término no es nulo para la dispersión entre valles. La energía final (y la estructura de bandas) producen directamente el módulo del nuevo momento k'. En este punto, solo es necesario elegir una nueva dirección (o ángulo) para la partícula dispersada. En algunos casos simples como la dispersión de fonones y una relación de dispersión parabólica, el ángulo de dispersión es aleatorio y se distribuye uniformemente en la esfera de radio k'. Usando coordenadas esféricas, el proceso de elección del ángulo es equivalente a elegir aleatoriamente dos ángulos y . Si el ángulo se distribuye con una distribución , entonces para una distribución uniforme de ángulos, la probabilidad de elegir un punto de la esfera es
Es posible, en este caso, separar las dos variables. Integrando sobre y luego sobre , se obtiene
Los dos ángulos esféricos pueden entonces elegirse, en el caso uniforme, generando dos números aleatorios 0 < r 1 , r 2 < 1 tales que
La tendencia actual de reducir el tamaño de los dispositivos semiconductores ha obligado a los físicos a incorporar cuestiones de mecánica cuántica para adquirir un conocimiento profundo del comportamiento de los dispositivos. Simular el comportamiento de dispositivos a escala nanométrica requiere el uso de un modelo de transporte cuántico completo , especialmente para los casos en los que no se pueden ignorar los efectos cuánticos. Sin embargo, esta complicación se puede evitar en el caso de dispositivos prácticos como el MOSFET moderno , empleando correcciones cuánticas dentro de un marco semiclásico. El modelo de Monte Carlo semiclásico se puede emplear entonces para simular las características del dispositivo. Las correcciones cuánticas se pueden incorporar a un simulador de Monte Carlo simplemente introduciendo un término de potencial cuántico que se superpone al potencial electrostático clásico visto por las partículas simuladas. La figura de al lado representa gráficamente las características esenciales de esta técnica. Los diversos enfoques cuánticos disponibles para la implementación se describen en las siguientes subsecciones.
La ecuación de transporte de Wigner constituye la base para la corrección cuántica basada en Wigner. [ cita requerida ]
donde k es el momento del cristal, V es el potencial clásico, el término en el lado derecho es el efecto de la colisión, el cuarto término en el lado izquierdo representa efectos mecánicos cuánticos no locales. La ecuación de transporte de Boltzmann estándar se obtiene cuando los términos no locales en el lado izquierdo desaparecen en el límite de variaciones espaciales lentas. La ecuación de transporte de Boltzmann simplificada (para ) corregida cuánticamente se convierte entonces en
donde el potencial cuántico está contenido en el término (debe ser un error: nunca se mencionó).
Este método de corrección cuántica fue desarrollado por Feynman y Hibbs en 1965. [ cita requerida ] En este método, el potencial efectivo se deriva calculando la contribución a la integral de trayectoria de las fluctuaciones cuánticas de una partícula alrededor de su trayectoria clásica. Este cálculo se lleva a cabo mediante un método variacional utilizando un potencial de prueba de primer orden. El potencial clásico efectivo en el punto promedio de cada trayectoria se convierte entonces en
Este método consiste en resolver periódicamente una ecuación de Schrödinger en una simulación con la entrada como potencial electrostático autoconsistente. Los niveles de energía exactos y las funciones de onda relacionadas con la solución del potencial electrostático se utilizan para calcular el potencial cuántico. La corrección cuántica obtenida sobre la base de este método se puede visualizar mediante la siguiente ecuación
donde V schr es el potencial de corrección cuántica, z es la dirección perpendicular a la interfase, n q es la densidad cuántica de la ecuación de Schrödinger que es equivalente a la concentración de Monte Carlo convergente, V p es el potencial de la solución de Poisson, V 0 es el potencial de referencia arbitrario alejado de la región cuántica de modo que la corrección se vuelve nula en la región de comportamiento semiclásico. Aunque los potenciales mencionados anteriormente para la corrección cuántica difieren en su método de cálculo y sus supuestos básicos, cuando se trata de su inclusión en la simulación de Monte Carlo, todos se incorporan de la misma manera.