Résoudre l'équation de Schrödinger

Je me souviens de mon émerveillement lors de ma première rencontre avec l'équation de Schrödinger. C'était en math sup, au milieu des années 70. A l'époque, l'enseignement de la mécanique quantique en prépa était assez frustre, notre prof nous avait "balancé" cette équation, pratiquement au début de l'année : émerveillement mais subite angoisse ! Et pourtant j'adore cette équation, sa simplicité apparente dissimule un monde de perplexité !

Aujourd'hui dans cette page, je vais revenir sur le sujet et vous proposer une méthode numérique simple pour la résoudre. Et pour simplifier encore plus les choses, je me placerai dans un espace à une seule dimension spatiale. Les particules dont j'étudierai le mouvement sont non relativistes (sinon de toute manière, Schrödinger n'est pas applicable).

Voyons donc comment simuler le mouvement d'une particule en mécanique quantique. Pour illustrer le problème, considérons un électron placé dans une boite dite "quantique" car sa dimension est de l'ordre de la taille d'un atome, quelques nanomètres. Dans ces conditions pas question d'utiliser les lois de Newton : la dimension de la boite nous impose l'usage des règles de la mécanique quantique. Comment s'y prendre ?

La fonction d'onde

Les conséquences de la dualité onde-corpuscule

Dans la page de TangenteX consacrée à la dualité onde-corpuscule, j'ai abordé le comportement dual d'une particule, qui peut être décrit par un modèle corpusculaire ou bien ondulatoire. Rappelons-en quelques élements ici, car l'équation de Schrödinger est d'abord une équation d'onde...

En 1905, Einstein a proposé une hypothèse corpusculaire du rayonnement électromagnétique en étendant les résultats de Planck sur la quantification du rayonnement. Il postule que l'énergie est quantifiée selon la formule \( E = h \nu \) ou encore \( E = \hbar \omega \), en utilisant la fréquence ou la pulsation de l'onde électromagnétique.

D'autre part, à la même époque, Louis de Broglie sauve le modèle ondulatoire en postulant que chaque particule est accompagnée d'une onde (de nature indéterminée) et il établit une autre formule célébre : \( \hbar \vec{k} = \vec{p} \), \( \vec{k} \) étant le vecteur d'onde de ladite onde, et \(\vec{p} \) sa quantité de mouvement (ou impulsion). La longeur d'onde de cette onde est donnée par la relation \( \lambda = \dfrac{h}{p} = \dfrac{2 \pi}{k} \), qui relie le vecteur d'onde et la quantité de mouvement de la particule. Etonnant non, cette formule qui inclut deux variables de deux mondes différents, l'une ondulatoire et l'autre corpusculaire....

Donc de Broglie postule l'existence d'une onde harmonique qui se déplace dans la même direction que la particule associée, dont l'équation est  :

\( \Psi ( \vec{x},t) = A e^{i(\vec{x}.\vec{k} - \omega t)} \)

La fonction d'onde en mécanique quantique

Sur la base des idées de de Broglie, un concept fondamental de la mécanique quantique fut créé : la fonction d'onde, encore nommée fonction d'état. Cette fonction décrit entièrement l'état d'une particule quantique. Dans un espace unidimensionnel, on la note \( \Psi(x,t) \) et elle est définie par :

\( \Psi(x,t) = \Psi_0 e^{i(kx - \omega t)} \)

C'est une fonction à variables complexes. J'insiste sur un point capital : il s'agit d'un objet mathématique, qui n'a aucun sens physique ! Dans l'interprétation de Copenhague, on postule que cette fonction décrit une amplitude de probabilité, qui permettra de calculer la densité de probabilité de détection d'une particule à un endroit donné en un lieu donné. A ce propos, beaucoup d'ouvrages utilisent le terme "densité de probabilité de présence" : je n'aime pas trop car la notion de "présence" d'une particule en un lieu et instant donnés n'a pas de sens si elle n'est pas détectée, au sens "mesure" du terme (voir la notion de réduction du paquet d'ondes, ici par exemple).

La densité de probabilité de détection se calcule simplement comme le produit de la fonction d'onde par son conjugué \( \Psi \Psi^* \), qui comme vous le savez est aussi égal au carré de son module \( ||\Psi||^2 =  \Psi \Psi^* \). Dans l'espace monodimensionnel, l'unité de la densité de probabilité de détection est le m-1.

Quelques formules utiles pour les calculs à propos des relations entre le vecteur d'onde k et la longueur d'onde associée : \( k = \dfrac{2 \pi}{\lambda}\) et \( \omega = \dfrac{\hbar k^2}{2m_e}\).

Les conséquences physiques de la notion de fonction d'onde et de densité de probabilité de détection sont immenses, comme par exemple la disparition de la notion de trajectoire d'une particule qui n'a aucun sens en mécanique quantique. Mais dans cette page, je me limite à la résolution numérique de l'équation de Schrödinger, je ne les aborderai donc pas.

Puisque nous parlons de probabilités, nous pouvons en déduire que la somme des probabilités de détection doit être égale à 1, c'est à dire que si une particule est prisonnière dans notre boite quantique unidimensionnelle d'extension comprise entre 0 et L, alors nous devons avoir :

\( \int_0^L  ||\Psi||^2 dx = 1 \)

ce qu'on appelle la condition de normalisation de la fonction d'onde. Cela impose que la fonction d'onde \( \Psi(x,t) \) soit dite de "carré sommable", autrement dit que l'intégrale ci dessus ne doit pas diverger pour avoir un sens physique. Les physiciens ont horreur des intégrales qui divergent...

La représentation d'un électron par un paquet d'ondes

La fonction d'onde tel qu'elle est définie par de Broglie est une onde plane progressive harmonique, comme vous avez pu en rencontrer en optique ou en électromagnétisme. Si l'on veut lui donner une signification physique en mécanique quantique, on se heurte vite à un problème : son module est le même partout ! Ce qui est gênant, car cela signifie que la densité de probabilité de détection de notre électron est la même partout. Et notre électron ne peut pas être partout avec la même probabilité !

Pour résoudre ce problème, on améliore la fonction d'onde définie par de Broglie, en ne considérant plus seulement une onde plane, mais une superposition d'ondes planes monochromatiques, un paquet d'ondes. Après tout, l'équation de Schrödinger est linéaire et donc une combinaison linéaire de solutions est aussi une solution.

On construit donc un paquet d'ondes en additionnant des ondes monochromatiques de fréquences différentes mais proches, de telle sorte que notre fonction d'onde soit de carré sommable. La fonction \( \Psi(x,t) \) devient donc, dans le domaine borné de tout à l'heure :

\( \Psi(x,t) = \int_0^L A(k) e^{i(kx - \omega t)} dk \)

la fonction A(k) définissant le spectre de fréquences de notre paquet d'ondes.

On choisit une enveloppe gaussienne pour distribuer les fréquences dans notre paquet d'ondes car nous obtenons alors une densité de probabilité de détection de forme gaussienne, pour laquelle la valeur la plus probable est la valeur moyenne, ce qui simplifie la vie. Il y a d'autres raisons, mais je ne vais pas m'apesentir sur le sujet.

Les paquets d'ondes possèdent des propriétés très intéressantes, sur lesquelle je ne m'étendrai pas ici. Signalons simplement (sic !) que l'on en déduit les inégalités d'Heisenberg.... Vous pourrez observer certaines de ces propriétés en manipulant, par exemple l'étalement dans le temps du paquet d'ondes.

Rappelons l'expression de l'énergie moyenne transportée par un paquet d'ondes : \( E = \dfrac{\hbar^2 k^2}{2 m_e} \) ou encore exprimée en fonction de la longueur d'onde \( E = \left({\dfrac{h}{\lambda}} \right)^2  \dfrac{1}{2m_e}\).

L'équation de Schrödinger "time-dependent"

Soit une particule, un électron pour fixer les idées, qui se déplace dans un espace unidimensionnel. Pour étudier son mouvement, nous allons fixer un référentiel. L'apparition de ce mot "référentiel" devrait vous interpeler : la relativité n'est pas loin ! On ne va pas se compliquer la vie et on considérera que la vitesse de la particule est petite devant c, un électron non relativiste. Pour rappel, si la particule est relativiste, ne surtout pas utiliser Schrödinger...

Cet électron bouge dans un espace où il peut subir des interactions avec d'autres systèmes physiques. Dans ce cas, on modélisera ces interactions par une énergie potentielle notée U(x,t). Si l'électron ne subit aucune interaction, U(x,t) sera nulle et notre électron sera dit "libre".

L'équation, analogue à notre deuxième loi de Newton, qui décrit le mouvement de l'électron dans l'espace et le temps est l'équation de Schrödinger. Plus précisément, elle décrit l'évolution de la fonction d'onde qui elle-même décrit l'état quantique de l'électron. Sa forme unidimensionnelle s'écrit :

\( i\hbar \dfrac{\partial \Psi (x,t)}{\partial t} = -\dfrac{\hbar ^2}{2m_e} \dfrac{\partial^2 \Psi (x,t)}{\partial x^2} + U(x)\Psi(x,t) \)

avec me la masse de l'électron et \( \hbar \) la constante de Planck réduite. Je choisis d'utiliser une fonction potentiel U stationnaire (indépendante du temps).

C'est une équation linéaire de premier ordre par rapport au temps. Linéaire, c'est à dire que toute combinaison linéaire de solutions particulières de l'équation est aussi solution de l'équation. De premier ordre par rapport au temps, c'est à dire que la connaissance de l'état quantique de la particule à un instant donné, considéré comme l'instant initial t0, permet de calculer son état quantique à n'importe quel instant postérieur à t0.

Le problème

L'énoncé

Revenons à notre problème : on enferme un électron dans une boite quantique, c'est à dire qu'on se cantonne dans un espace fermé par des barrières de potentiel infini. Nous resterons dans un espace unidimensionnel, ce qui nous donne le schéma suivant de boîte :

Boite quantique

Dans cette boîte, nous déposons un électron et nous voulons déterminer son mouvement. Pour l'instant, notre électron ne rencontre aucune barrière de potentiel dans sa boîte. Il est libre... dans sa boîte ! La fonction U(x) est donc constante et nulle dans la boîte.

Pour ce faire, nous allon mobiliser 4 outils : la notion de paquet d'ondes, l'équation de Schrödinger, l'algorithme FDTD (que nous avons déjà rencontré ici) et Python.

La simulation de l'évolution d'un électron libre U(x) = 0

Créer le paquet d'ondes d'enveloppe gaussienne

Simuler un électron avec un paquet d'ondes

Nous allons donc simuler l'électron par une fonction d'onde d'enveloppe gaussienne. Pour faciliter les calculs, j'ai décomposé, comme c'est l'habitude, la fonctions d'onde qui est une fonction complexe, en sa partie réelle et sa partie imaginaire, que je vais calculer séparément. Ainsi l'équation :

\(\Psi = e^{-\dfrac{1}{2}{\left(\dfrac{(x-x_0)}{\sigma}\right)}^2}e^{i\dfrac{2\pi(x-x_0)}{\lambda}} \)

qui définit le paquet d'ondes à enveloppe gaussienne, avec \( \sigma \) la largeur du paquet et x0 l'abscisse initiale du centre du paquet d'ondes, devient en passant l'exponentielle sous forme trigonométrique :

\( \Psi_{Reel} =  e^{-\dfrac{1}{2}{\left(\dfrac{(x-x_0)}{\sigma}\right)}^2} cos \left(\dfrac{2\pi(x-x_0)}{\lambda}\right) \)

\( \Psi_{Imag} =  e^{-\dfrac{1}{2}{\left(\dfrac{(x-x_0)}{\sigma}\right)}^2} sin\left(\dfrac{2\pi(x-x_0)}{\lambda}\right) \)

Il restera à normaliser le paquet initial pour satisfaire la condition :

\( \int_0^L \left( \Psi_{Reel}^2 +  \Psi_{Imag}^2  \right)dx = 1 \).

Dans ma simulation j'ai choisi les valeurs suivantes pour les paramètres physiques :

Vous pourrez faire varier la longueur d'onde associée à l'électron, et donc son énergie, en observant ce qui se passe. Vous pouvez aussi modifier la largeur du paquet, en faisant les mêmes observations.

Le code Python

Il est contenu dans le script SchrodingePaquetOndes.py. Ce script commence classiquement par l'import des packages utiles et par la déclaration des différents paramètres de la simulation. J'utilise ici le package scipy.constants qui contient la plupart des constantes physiques, ce qui évite d'aller les chercher dans les bouquins.

Je commence par calculer la valeur de la partie réelle et de la partie imaginaire de la fonction d'onde initiale:
Psi_Real = exp(-0.5*((x-x0)/sigma)**2)*cos(DeuxPi*(x-x0)/Lambda)
Psi_Imag = exp(-0.5*((x-x0)/sigma)**2)*sin(DeuxPi*(x-x0)/Lambda)

puis je procède à sa normalisation et au calcul de la densité de probabilité de présence \( ||\Psi||^2 \), comme décrit plus haut :

PsiPsiC = Psi_Real**2 + Psi_Imag**2
C = simps(PsiPsiC,x)
Psi_Real = Psi_Real/sqrt(C)
Psi_Imag = Psi_Imag/sqrt(C)
Psi_Prob = Psi_Real**2 + Psi_Imag**2

J'ai employé ici la méthode de Simpson, tirée du package integrate, pour calculer l'intégrale sous le paquet d'ondes sur le domaine [0,L]. J'aurais aussi pu employer la méthode des trapèzes (moins précise) ou Romberg (plus couteuse).

Puis j'ai calculé l'énergie transportée par le paquet, c'est à l'énergie de l'électron simulé et j'indique l'énergie théorique de cet électron :
En = zeros(Nx)
En[1:-1] = a1*(Psi_Real[1:-1] - 1j*Psi_Imag[1:-1]) \
           *(Psi_Real[2:] - 2*Psi_Real[1:-1] + Psi_Real[:-2] \
           + 1j*(Psi_Imag[2:] - 2*Psi_Imag[1:-1] + Psi_Imag[:-2]))                    
Ebarre = simps(En,x)
print "Energie theorique : " + str(Ec) + " eV"
print "Energie moyenne du paquet : " + str(Ebarre.real) + " eV"

L'énergie théorique calculée est de 66,84 eV et l'énergie moyenne du paquet est de 67,07 eV.

Visualisation d'un paquet d'ondes initial

Les instructions d'affichage n'ont rien d'original !

Les résultats

Sur le schéma suivant, on observe la fonction d'onde d'un paquet d'ondes gaussien (en bleu la partie réelle et en rouge la partie imaginaire) :

Fonction d'onde gaussienne

et ici la courbe de densité de probabilité de détection :

Densité de probabilité de détection

On constate qu'il s'agit d'une distribution gaussienne centrée sur l'abscisse d'origine du paquet d'ondes.

Résoudre Schrödinger avec l'algorihme FDTD 1D

Sa résolution analytique

L'équation de Schrödinger 1D avec U(x) = 0, que l'on écrit \( i\hbar \dfrac{\partial \Psi (x,t)}{\partial t} = -\dfrac{\hbar ^2}{2m} \dfrac{\partial^2 \Psi (x,t)}{\partial x^2} \) est parfaitement intégrable analytiquement ! D'ailleurs, je vous invite à la comparer avec une autre équation que nous avons déjà rencontré sur TangenteX : \( \dfrac{\partial T}{\partial t} = D \dfrac{\partial^2 T}{\partial x^2} \), l'équation de diffusion thermique ! Dans les deux cas, la transformée de Fourier est l'outil privilégié pour ces calculs.

Mais sur cette page, nous allons plutôt aborder le problème sous l'angle numérique, en utilisant l'algorithme FDTD, sachant que d'autres méthodes sont possibles.

La discrétisation de Schrödinger selon la méthode FDTD 1D

Nous avons déjà rencontré l'algorithme FDTD (Finite Differences Time Domain) 1D lorsque je vous ai proposé une simulation de la propagation d'une onde électromagnétique.

La première étape consiste à discrétiser l'équation de Schrödinger. Reprenons donc notre équation de Schrödinger time-dependent dans un espace unidimensionnel :

\( i\hbar \dfrac{\partial \Psi (x,t)}{\partial t} = -\dfrac{\hbar ^2}{2m} \dfrac{\partial^2 \Psi (x,t)}{\partial x^2} + U(x)\Psi(x,t) \)

La fonction \( \Psi(x,t) \) est une fonction complexe. Je vais la décomposer en sa partie réelle \( \Psi_R(x,t)\) et sa partie imaginaire \( \Psi_I(x,t) \) et écrire deux équations l'une portant sur sa partie réelle et l'autre sur sa partie imaginaire, ce qui, en réorganisant notre équation de Schrödinger pour isoler à gauche la dérivée partielle temporelle, nous donne le couple d'équations :

\( \dfrac{\partial \Psi_R (x,t)}{\partial t} = -\dfrac{\hbar}{2m} \dfrac{\partial^2 \Psi_I (x,t)}{\partial x^2} + \dfrac{1}{\hbar}U(x)\Psi_I(x,t) \)

\( \dfrac{\partial \Psi_I (x,t)}{\partial t} = \dfrac{\hbar}{2m} \dfrac{\partial^2 \Psi_R (x,t)}{\partial x^2} - \dfrac{1}{\hbar}U(x)\Psi_R(x,t) \)

Comme d'habitude je vais discrétiser les dérivées partielles en faisant un DL d'ordre 1 :

\( \dfrac{\partial \Psi(x,t)}{\partial t} \approx \dfrac{\Psi(x, t+\Delta t) - \Psi(x,t)}{\Delta t} \)

\( \dfrac{\partial^2 \Psi(x,t)}{\partial x^2} \approx  \dfrac{\Psi(x+\Delta x,t) - 2\Psi(x,t) + \Psi(x - \Delta x,t)}{{\Delta x}^2}\)

Reste à reporter ces valeurs dans les équations et à appliquer l'algorithme de Yee pour obtenir les équations discrètes qui vont nous servir à coder le FDTD 1D. Je vous renvoie à la page sur la propagation des OEM pour plus de détails.

Nous obtenons finalement les équations discrètes, avec i l'indice sur le vecteur x :

\( \Psi_R(i)  = \Psi_R(i) - \dfrac{\hbar \Delta t }{2m{\Delta x}^2}\left( \Psi_I(i+1) - 2\Psi_I(i) + \Psi_I(i-1)  \right)  + \dfrac{e \Delta t}{\hbar}  U(i) \Psi_I(i) \)

\( \Psi_I(i)  = \Psi_I(i) + \dfrac{\hbar \Delta t }{2m{\Delta x}^2}\left( \Psi_R(i+1) - 2\Psi_R(i) + \Psi_R(i-1)  \right)  - \dfrac{e \Delta t}{\hbar}  U(i) \Psi_R(i) \)

Une petite remarque : j'ai introduit la charge élémentaire e dans la valeur du coefficient affectant l'énergie potentielle U : c'est parce que j'exprime dans mon code U en eV, il faut donc effectuer la conversion eV/Joule.

Attention aux conditions de stabilité numérique de l'algorithme de Yee, qui contraignent le pas spatial \( \Delta x \) et temporel \( \Delta t \). Dans mon code, j'ai défini \( a_2 = \dfrac{\hbar \Delta t }{2m{\Delta x}^2} = 0.1 \), il faut bien faire un choix, ce qui me donne la contrainte sur \( \Delta t \) suivante, codée dans la variable dt : dt = a2*2*m_e*dx**2/hbar.

Le code Python

Il est contenu dans le script SchrodingeFDTD_1D.py. Le début du script est une copie du script SchrodingePaquetOndes.py. Vous pourrez là aussi faire varier la longueur d'onde associée à l'électron, et donc son énergie, en observant ce qui se passe. Vous pourrez aussi modifier la largeur du paquet, en faisant les mêmes observations. Les valeurs de dx et dt sont fixées pour assurer la stabilité numérique de l'algorithme. Les conditions aux limites \( \Psi(0,t) = 0 \) et \( \Psi(L,t) = 0 \) sont fixées lors de l'initialisation des buffers de calcul Psi_Real, Psi_Imag et Psi_Prob.

Dans ce premier code, je fixe la fonction Potentiel U à zero avec l'instruction U = zeros(Nx).

Bien que l'algorithme FDTD soit d'apparence rébarbative, son codage est très simple, surtout en Python. Au lieu d'utiliser des boucles, j'exploite le slicing de Python qui permet d'améliorer le temps de calcul de plusieurs ordres de grandeur. Ce n'est pas immédiat, mais avec un peu de réflexion, c'est très efficace.

De fait, le code tient en trois lignes, y compris le calcul de la densité de probabilité de détection courante:

Psi_Real[1:-1] = Psi_Real[1:-1] - a2*(Psi_Imag[2:] - 2*Psi_Imag[1:-1] + Psi_Imag[:-2]) + a3*U[1:-1]*Psi_Imag[1:-1]
Psi_Imag[1:-1] = Psi_Imag[1:-1] + a2*(Psi_Real[2:] - 2*Psi_Real[1:-1] + Psi_Real[:-2]) - a3*U[1:-1]*Psi_Real[1:-1]
Psi_Prob[1:-1] = Psi_Real[1:-1]**2 + Psi_Imag[1:-1]**2

Visualiser l'évolution de l'électron

J'utilise ici les fonctions standard de pyplot. Il est possible d'utiliser d'autres outils comme le package matplotlib.animation, qui permet d'animer des courbes, ou Visual Python, mais j'ai préféré rester simple. Un détail pratique : cela ne marche pas dans l'environnement Anaconda/Spyder. J'utilise donc utilisé VIDE.

Le code est très simple :

    if t % PasAff == 0:       
        linePsiR.set_ydata(Psi_Real)
        linePsiI.set_ydata(Psi_Imag)
        linePsiP.set_ydata(Psi_Prob/1.e5)
        pause(0.1)

Les objets linePsiR, linePsiI et linePsiP sont initialisés lors du tracé du paquet d'ondes initial. Les courbes sont modifiées en utilisant la méthode set_ydata(). Une méthode à ne pas oublier lorsqu'on veut animer des courbes simplement avec pyplot.

Le test if t % PasAff == 0 permet de n'afficher qu'une courbe toutes les PasAff pas de calcul, pour plus de fluidité. La valeur de la pause (en secondes) est ajustable en fonction de votre besoin.

Les résultats

Voici une vue à t = 2 000 de l'animation. La courbe bleue représente la partie réelle de la fonction d'onde, en rouge sa partie imaginaire. La courbe verte représente la densité de probabilité de détection de l'électron (\( \Psi \Psi^* \))  :

Schrodinger FDTD 1D V(x,t) =0

Si vous suivez l'évolution du paquet d'ondes avant sa réflexion sur les parois de la boite quantique, vous observerez plusieurs choses :

Voyons ces différents points.

La propagation du paquet d'ondes

Le paquet d'ondes se propage ! Nous pouvons donc introduire la notion de vitesse du paquet, une superposition d'ondes. C'est ce que les physiciens appellent une vitesse de groupe. La vitesse de groupe du paquet est égale à \(v_g = \dfrac{h}{m_e\lambda} \), ce qui est vérifiable graphiquement en relevant la position du maximum de la courbe de densité de probabilité de détection en fonction du temps. Vous aussi pourrez vérifier que plus la longueur d'onde moyenne du paquet \( \lambda \) est petite, i.e. plus son énergie est grande, et plus le paquet se déplace vite.

Nous observons également que lorsque le paquet atteint la barrière de potentiel infini qui délimite la boite quantique, il se réfléchit, comme le ferait un paquet d'ondes électromagnétique.

L'étalement du paquet

A t=0 le paquet possède une quantité de mouvement (on emploie fréquement le terme impulsion) fixe, et qui restera constante pendant toute l'évolution du paquet.

Les relations d'inégalité d'Heisenberg nous permettent d'écrire que \( \Delta x \Delta p_x \ \ge \dfrac{\hbar}{2}\). Or, l'impulsion est constante et donc c'est l'extension spatiale du paquet qui augmente : il s'étale.

La variation d'amplitude du paquet

Cette variation découle de l'étalement du paquet : son énergie demeure constante et son extension spatiale augmente, donc l'amplitude du paquet diminue. Rappelez-vous, on calcule l'énergie moyenne du paquet par l'intégration sous l'enveloppe du paquet...

Simulation de l'évolution de l'électron avec U(x) non nulle

Cas d'une marche de potentiel

Essayons maintenant de suivre l'évolution de notre électron si nous disposons dans sa boite une marche de potentiel, dont nous pourrons régler la hauteur. En fait, rien n'est plus simple : aucun changement d'algorithme ni de code à faire dans la partie concernant la résolution de l'équation de Schrödinger. Il suffit de définir une nouvelle fonction U puis de l'afficher sur le graphe, en faisant attention aux échelles.

Définir la marche de potentiel

Le code est très simple. Je fixe la hauteur de la marche :

U0 = 80                # en eV

puis je construit un buffer de zéros, dont je modifie la valeur à partir de la position de la marche que j'ai fixé à la moitié de la boite :

U = zeros(Nx)
U[Nx/2:] = U0          # définition d'une marche de potentiel

Voilà, c'est tout.... Pour modifier la hauteur de la marche, il suffit de modifier la valeur de U0, en tenant compte toutefois de l'énergie de l'électron, soit approximativement 60 eV avec mes paramètres. Il est d'ailleurs intéressant de faire une marche moins haute que son énergie, ou beaucoup plus haute pour voir ce que cela donne. Il suffit de modifier la valeur U0 puis de relancer le script !

L'affichage de la marche

Comme j'ai envie d'afficher la marche sur le même graphe que la densité de probabilité de détection, il va falloir que je me préoccupe de la différence d'échelle et de grandeur, qui est considérable. J'utilise pour cela une fonction de matplotlib, twinx(), qui me permet de tracer la densité de probabilité avec son échelle, en indiquant son échelle sur l'axe des y à gauche, sur axe1, et la marche de potentiel sur l'axe des y de droite, avec son échelle, sur axe2. Pour afficher la marche :
lineEner, = axe2.plot(x*1.e9,U,'blue')

Les résultats

La hauteur de la marche de potentiel est fixée à 80 eV, alors que l'énergie moyenne de l'électron est de 67 eV. Lorsque l'électron percute la marche, voilà ce que devient la densité de probabilité de détection : elle est un peu chahutée. En fait, nous observons les interférences entre l'onde incidente, qui traverse la marche et l'onde qui se réfléchit sur la marche :

Schrodinger FDTD 1D V(x,t) = marche

Mais le plus étonnant, c'est ce qui se passe après la percussion sur la marche :

Schrodinger FDTD 1D V(x,t) = marche

En mécanique classique, il est impossible que l'électron franchisse la marche de potentiel, son énergie est trop faible pour cela. Et pourtant, nous constatons que la densité de probabilité de détection de l'électron à l'intérieur de la marche n'est pas nulle ! C'est cette caractéristique qui permet d'expliquer certains phénomènes comme la radioactivité alpha, phénomènes absolument incompréhensible en physique classique.

Nous observons aussi que :

Cas d'une barrière de potentiel

Définir la barrière de potentiel

Le principe de définition est le même que pour une marche : on déclare un buffer U à zéros, puis on le modifie en définissant la barrière. Ici, il nous faut deux paramètres : U0 sa hauteur et EppBar son épaisseur. Ce qui nous donne le code

EppBar = 100                    # largeur de la barrière en nombre de pas dx
U[Nx/2:Nx/2+EppBar] = U0        # définition d'une barrière de potentiel

L'affichage de la barrière

Comme pour la marche, on ne change rien !

Les résultats

La hauteur de la barrière de potentiel est fixée à 80 eV, alors que l'énergie moyenne de l'électron est de 67 eV. Son épaisseur est fixée à 0,15 nanomètres. Lorsque l'électron percute la barrière, voilà ce que devient la densité de probabilité de détection :

Schrodinger FDTD 1D V(x,t) = barrière

Mais le plus étonnant, c'est ce qui se passe au delà de la barrière :

Schrodinger FDTD 1D V(x,t) = barrière

La densité de probabilité de détection de l'électron derrière la barrière de potentiel n'est pas nulle. Cela signifie que la probabilité que notre électron franchisse la barrière de potentiel n'est pas nulle, même si son énergie est plus petite que la hauteur de la barrière. C'est le fameux effet tunnel, que l'on retrouvera aussi en étudiant la radioactivité alpha.

Essayez de faire varier l'énergie de l'électron, en diminuant sa longeur d'onde (Lambda), ou bien faites varier la largeur (EppBar) et la hauteur de la barrière de potentiel (U0), et observez les résultats.

Les scripts Python

Les scripts Python étudiés dans cette page sont disponibles dans le package Schrodinger1D.zip :

Contenu et design par Dominique Lefebvre - www.tangenteX.com septembre 2018 Licence Creative Commons Contact :

Cette œuvre est mise à disposition selon les termes de la Licence Creative Commons Attribution - Pas d’Utilisation Commerciale - Pas de Modification 3.0 France.