Oscillateur harmonique en une dimension

UPMC, Paris. (P Reinhardt, LCT)

Clear["Global`*"]

Dans un espace à une dimension, on considère un système formé de deux particules de masses (Subscript[m, 1]) et (Subscript[m, 2]) évoluant dans un potentiel harmonique V(x)=(1/2)k (x^2). Cela revient à étudier le mouvement d'une seule particule fictive de masse réduite μ=(Subscript[m, 1])(Subscript[m, 2])/((Subscript[m, 1])+(Subscript[m, 2]))  soumise au potentiel V(x).

A)  En mécanique classique :

Rappeler l'équation du mouvement de la particule et résoudre cette équation différentielle avec Mathematica avec les conditions initiales x(t=0)=(Subscript[x, 0]) et x'(t=0)=0.

solution=DSolve[{μ x''[t]+k x[t]==0,x[0]==x0,x'[0]==0},x[t],t]

Tracer le résultat (en fonction de  t ) sur l'intervalle de temps [-10,10] pour (Subscript[x, 0])=1, μ=0.5  pour différentes valeurs de k.

Sol[t_,k_]=x[t]/.solution[[1]]

x0=1;μ=0.5;

Plot[{Sol[t,1],Sol[t,5],Sol[t,10],},{t,-10,10}]

Pour obtenir une densité de probabilité de présence en fonction de x il faut intégrer le mouvement sur une demi-période :

  T=(∫_x0^(-x0))(dt/dx)dx= -(∫_(-x0)^x0)(dt/dx)dx
 Avec la substitution x(t) = x0 cos(ωt) on trouve avec  x0 sin(α)=x0(Sqrt[1-(cos^2)(α)]) pour la derivée dx/dt=-ω (Sqrt[(x0^2)-(x^2)(t)]), et l'integrale devient 
   T=(1/ω)(∫_(-x0)^x0)(1/(Sqrt[(x0^2)-(x^2)]))dx =(π/ω)     
  La densité de probabilité de présence classique pour trouver la particule à l'endroit  x est alors 
 
 ProbClass[x_] := (1/π) / Sqrt[x0^2 - x^2]

Intégrer cette densité de probabilité de présence entre -x0 et x0. Tracer la courbe avec x0 = 2.

ProbClass[x_]:=(1/π)/Sqrt[x0^2-x^2]

Integrate[ProbClass[x],{x,-x0,x0}]

x0=2;Plot[ProbClass[x],{x,-x0,x0},PlotRange->{{-x0,x0},{0,1}}]

B) En mécanique quantique :  
L'équation de Schrodinger pour l'oscillateur harmonique s'écrit

-((ℏ^2)/(2μ))(((d^2)ψ)/(d (x^2)))+(1/2)μ (ω^2)(x^2)ψ=E ψ
avec ω=(Sqrt[k/μ]).On peut effectuer le changement de variable suivant :
ξ=(Sqrt[((μω)/ℏ)])x      et    ϵ = (2/(ℏω))E
qui conduit à l'équation suivante :

-(((d^2)ψ)/(d (ξ^2)))+(ξ^2)ψ=ϵ ψ

Si  ξ->∞ le terme en (ξ^2) devient prépondérant et le terme en ϵ devient négligeable. Cela nous permet de dire que l'on peut réécrire ψ de la façon suivante:

ψ(ξ)=ϕ(ξ)(e^(-((ξ^2)/2)))

 Si on remplace ψ(ξ) dans l'équation différentielle précédente on obtient finalement l'équation différentielle suivante pour ϕ(ξ) :
 
(((d^2)ϕ(ξ))/(d (ξ^2)))-2ξ ((d ϕ(ξ))/(d ξ))+(ϵ-1) ϕ(ξ)=0

Résoudre cette derniere équation (par calcul formel)

DSolve[Phi''[xi]-2 xi Phi'[xi]+(e-1)Phi[xi]==0,Phi[xi],xi]

Vérifiez dans l'aide les fonctions obtenus. Que pouvez-vous en déduire sur le premier terme de la fonction HermiteH ? Sur le terme C[2] ?

Le premier argument de HermiteH doit être un entier, ce qui fixe (ϵ - 1)/2 =  n avec n = 0,1, ... .  De plus, ψ(ξ) ne doit pas diverger quand ξ->∞ ;  or la fonction Hypergeometric1F1 conduit à une divergence, donc  le deuxième terme doit disparaitre : C[2] = 0

Donner l'expression de la fonction d'onde ψ(ξ) complète en fonction de n. (Définition d'une fonction Ψ(ξ,n) pour désigner (Subscript[Ψ, n])(ξ) )

Psi[xi_,n_]:=HermiteH[n,xi] Exp[-(xi^2)/2]

Psi[x,n]

Calculez la norme de ces fonctions ( (Sqrt[(∫_(-∞)^(+∞))(SuperStar[ψ])(ξ)ψ(ξ)dξ])    ) en fonction de n.

N0=1/Sqrt[Integrate[Psi[xi,0]^2,{xi,-∞,∞}]]
N1=1/Sqrt[Integrate[Psi[xi,1]^2,{xi,-∞,∞}]]
N2=1/Sqrt[Integrate[Psi[xi,2]^2,{xi,-∞,∞}]]
N3=1/Sqrt[Integrate[Psi[xi,3]^2,{xi,-∞,∞}]]
N4=1/Sqrt[Integrate[Psi[xi,4]^2,{xi,-∞,∞}]]
N5=1/Sqrt[Integrate[Psi[xi,5]^2,{xi,-∞,∞}]]

Affichez les premières 5 fonctions corrigées de leur coefficient de normation sous leur forme la plus simple. Puis, pour chacune d'entre elles, tracez la fonction et son carré (probabilité de présence) sur  l'intervalle [-6 ,6]

Psi0[x_]:=N0 Psi[x,0]
Psi1[x_]:=N1 Psi[x,1]
Psi2[x_]:=N2 Psi[x,2]
Psi3[x_]:=N3 Psi[x,3]
Psi4[x_]:=N4 Psi[x,4]
Psi5[x_]:=N5 Psi[x,5]

Plot[{Psi0[x]+1/2,Psi1[x]+3/2,Psi2[x]+5/2,Psi3[x]+7/2,Psi4[x]+9/2,Psi5[x]+11/2,1/2 x^2},{x,-5,5}]

Tracez la probabilité de présence pour n = 120 et comparez au résultat classique. Mettez (Subscript[x, 0])=(Sqrt[2n+1]).

n=120
x0=Sqrt[2n+1]
ProbClass[x_]:=(1/π)/Sqrt[x0^2-x^2]
NQuant=1/Integrate[Psi[x,n]^2,{x,-∞,∞}]
ProbQuant[x_]:=NQuant*Psi[x,n]^2
Plot[{ProbClass[x],ProbQuant[x]},{x,-1.3x0,1.3x0},PlotRange->{{-1.3 x0,1.3x0},{0,3.5/x0}}]