Equation de schrodinger : calcul des états propres dans un puits (méthode type euler)
Description
Ce code code permet de calculer les etats propres d'une particule quantique (quanton) dans un puit de potentiel (symetrique et "infini").
Il se decompose en 4 etapes principales.
_Premiere etape: il vous demande d'entrer n (numero du mode que vous desirez) au clavier
_Deuxieme etape: il definit la parite de la fonction d'onde correspondante
_Troisieme etape: apres avoir fait diverses operations mineures il calcule le fondamental et tous les modes inferieurs au n que vous avez entrer (mais il ne les sauvegarde pas!), et affiiche l'energie correspondant a chaque mode.
_Quatrieme etape: il inscrit ensuite les valeurs de la fonction d'onde que vous avez demander dans un ficher nomme Energie.plt
Voila il ne vous reste plus qu'a admirer le resultal et tracant la fonction avec un logiciel comme gnuplot par exemple.
Le potentiel par defaut est un potentiel en x². Vous pouvez essayer diverses formes (puit carre ou bizaroides) ils marchent tous du moment qu'ils sont suffisament abrupts, profonds et symetriques (pour le changer faites le dans le source...).
Pour plus de detail sur son fonctionnement lisez les commantaires, j'ai essayer d'etre le plus clair possible.
Et si vous avez un doute faites un peu de theorie, ca n'a jamais tuer personne... Pour ceux que ca interressent les livres de LE BELLAC ("Physique quantique") et LEVY-LEBLOND ("Quantique") sont assez exaustifs pour une premiere approche (mais allez y doucement quand meme ;) )
Nb: La mise en page est assez moche si vous essayer de lire le texte avec un editeur classique, je vous conseille donc de le regarder sous vi (pour les puristes) ou xemacs (ou tout autre logiciel adapte a la programmation)
Je tiens a remercier specialement F.Hebert qui a coodonner ce projet.
J'essairai de donner une version plus complete de la facon dont marche ce programme ainsi que les ameliorations a y apporter ulterieurement.
Desole pour les maniaques de l'orthographe mais moi c'est pas mon truc alors....
Il se decompose en 4 etapes principales.
_Premiere etape: il vous demande d'entrer n (numero du mode que vous desirez) au clavier
_Deuxieme etape: il definit la parite de la fonction d'onde correspondante
_Troisieme etape: apres avoir fait diverses operations mineures il calcule le fondamental et tous les modes inferieurs au n que vous avez entrer (mais il ne les sauvegarde pas!), et affiiche l'energie correspondant a chaque mode.
_Quatrieme etape: il inscrit ensuite les valeurs de la fonction d'onde que vous avez demander dans un ficher nomme Energie.plt
Voila il ne vous reste plus qu'a admirer le resultal et tracant la fonction avec un logiciel comme gnuplot par exemple.
Le potentiel par defaut est un potentiel en x². Vous pouvez essayer diverses formes (puit carre ou bizaroides) ils marchent tous du moment qu'ils sont suffisament abrupts, profonds et symetriques (pour le changer faites le dans le source...).
Pour plus de detail sur son fonctionnement lisez les commantaires, j'ai essayer d'etre le plus clair possible.
Et si vous avez un doute faites un peu de theorie, ca n'a jamais tuer personne... Pour ceux que ca interressent les livres de LE BELLAC ("Physique quantique") et LEVY-LEBLOND ("Quantique") sont assez exaustifs pour une premiere approche (mais allez y doucement quand meme ;) )
Nb: La mise en page est assez moche si vous essayer de lire le texte avec un editeur classique, je vous conseille donc de le regarder sous vi (pour les puristes) ou xemacs (ou tout autre logiciel adapte a la programmation)
Source / Exemple :
/********************************************Inclusion des bibliotheques standard*******************************************/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
/*********************************************Definition des constantes generales******************************************/
#define C 1.0 //Constante de Planck Heu... ca serait cool de le prendre = a 1.... precision machine etc...
#define PI 3.141592654 //PI
#define l2 10 //Definition de la largeur du puit
#define h 0.001 //pas d'integration
/***********************************************Fonction definissant le potentiel*********************************************/
/*REMARQUE: plus le puit et profond plus je pourai resoudre de modes mais plus le programme sera lent.
Exemple une profondeur de -100 permet de resoudre "seulement" les 12 premier modes du puit carre*/
double V (double i)
{
return (-1000+(h*h*i*i/10));
}
main()
{
/******************************************Toutes les variables dont je me sert******************************************/
double *x; //Definition de la fonction d'onde
double *v; //Definition de la derivee premiere de la fonction d'onde
double r; //Definition du pas dans la recherche de l'energie
double E; //Definition de l'energie
double U; //Autre definition de l'energie utilisee pour la dichotomie
double e; //Autre definition de l'energie pour la recherche grossiere
double l; //Definition de la longueur de l'espace d'integration
int p,q,w; //Sert lors de la recherche des modes
double m; //Masse de la particule
int parite; //Definit la parite de la fonction
FILE *pf1,*pf2,*pf3; //Fichier dans lequel je vais stocker mes donnees
int i,j,k,n,b,lg,z; //Definition de(s) truc(s) a incrementer dans les boucles
/*******************************************Initialisation des variables**************************************************/
printf("Quel mode cherchez vous? (0 correspond au fondamental)\nn=");
scanf("%d",&n); //n definit quel mode l'on recherche
printf("\n\n");
b=n; //On sauvegarde n qqpart car on va l'utiliser et le changer apres
if (n%2==0)
parite=1; //On impose la parite de la fonction d'onde suivant le mode que l'on recherche
else
parite=0;
w=1;
for (i=0;i<n;i++)
{ //On calcule le nombre de fois que la fonction va changer
if (i%2==0) //de facon de diverger
w=w+2;
else
w=w+1;
}
l=20; //N sert a definir la longueur sur laquelle on integre
lg=(int)(l/h)+1; //Eje decoupe l'axe sur lequel j'integre en petits segments de longueur h
m=1; //Fmasse de la particule
E=((C*C)/(2*m))*((PI/(2*l2))*(PI/(2*l2)))*(n+1)*(n+1)-1000; //RDefinition de l'energie pour le puits carre "infini"
//ARemarque si on chage la profondeur du puit de potentiel il faut changer e aussi
e=-1000; //AEnergie minimale du puits
r=0.01; //LPas dans la recherche de l'energie
x=calloc(lg,sizeof(double)); //Initialisation des tableaux
v=calloc(lg,sizeof(double));
//si la fonction est paire on lui impose sa valeur en zero comme etant 1
//si ma fonction d'onde est paire la derivee est nulle en zero
//si la fonction est impaire on lui impose sa valeur en zero comme etant 0
//si la fonction est impaire sa derivee est 1
if (parite == 1)
{
x[0]=1.0;
v[0]=0.0;
}
else
{
x[0]=0.0;
v[0]=1.0;
}
pf1=fopen("Energie.plt","w"); //Fichiers dans lesquels je vais preserver mes resultats
n=0; //on initialise ces valeurs a zero a fin de faciliter
p=0; //l automatisation du programme
q=0;
for(i=0;i<w;i++)
{
/****************************Recherche primaire de l'energie du niveau fondamental********************************/
if (n%2==0)
parite=1; //On impose la parite de la fonction d'onde suivant le mode que l'on recherche
else
parite=0;
if (parite == 1)
{
x[0]=1.0;
v[0]=0.0; //Et l'on donne les conditions initialles correspondantes
}
else
{
x[0]=0.0;
v[0]=1.0;
}
if(i%2==0)
{
for(e;e<-98;e+=r)
{
for(z=0;z<lg-1;z++)
{
v[z+1]=v[z]+h*(-((2*m)/(C*C)))*(e-V((double)z))*x[z]; //Algorithme d'integration: sorte de runge-kutta
x[z+1]=x[z]+h*v[z+1]; //a l'ordre 1
}
if (x[(int)((l/h))] <0) //Si la facon de diverger change c'est que j'ai
{ //depasser l'energie exacte du fondamental
E=e; //Ici diverger est traduit par:
U=E-r; //"si la valeur de ma fonction au dela du potentiel
break; //vaut moins qu'une certaine valeur"
} //U sert a preserver une valeur pour laquelle
} //la divergence est differente afin d'obtenir
//les deux valeurs entre lesquelles je vais
//faire ma dichotomie
/************************Recherche dichotomique de l'energie du niveau fondamental********************************/
if (n%2==0)
parite=1; //On impose la parite de la fonction d'onde suivant le mode que l'on recherche
else
parite=0;
if (parite == 1)
{
x[0]=1.0;
v[0]=0.0; //Et l'on donne les conditions initialles correspondantes
}
else
{
x[0]=0.0;
v[0]=1.0;
}
for(j=0;j<100;j++)
{
for(k=0;k<lg-1;k++)
{
v[k+1]=v[k]+h*(-((2*m)/(C*C)))*(((E+U)/2)-V((double)k))*x[k]; //Je test la convergence pour une valeur
x[k+1]=x[k]+h*v[k+1]; //comprise entre E et U
}
if (x[(int)((l/h))] <0) //Si cette convergence ce fait de maniere identique a celle de E trouvée
//precedement alors c'est une energie plus proche de la valeur exacte
E=(E+U)/2;
if (x[(int)((l/h))]>0) //En revenche si la convergence ce fait de manire differente alors
//je peux augmenter ma borne inferieure (U) car mon energie n'est
U=(E+U)/2; //pas dans le demi intervalle inferieur de [E;U]
}
}
/*REMARQUE: la divergence tantot vers plus l'infini tantot vers moins l'infini explique la presence de la
deuxieme recherche primaire et dichotomique*/
if(i%2 != 0)
{
for(e;e<-98;e+=r)
{
for(z=0;z<lg-1;z++)
{
v[z+1]=v[z]+h*(-((2*m)/(C*C)))*(e-V((double)z))*x[z]; //Algorithme d'integration: sorte de runge-kutta
x[z+1]=x[z]+h*v[z+1]; //a l'ordre 1
}
if (x[(int)((l/h))] >0) //Si la facon de diverger change c'est que j'ai
{ //depasser l'energie exacte du fondamental
E=e; //Ici diverger est traduit par:
U=E-r; //"si la valeur de ma fonction au dela du potentiel
break; //vaut moins qu'une certaine valeur".
} //U sert a preserver une valeur pour laquelle
} //la divergence est differente afin d'obtenir
//les deux valeurs entre lesquelles je vais
//faire ma dichotomie
/************************Recherche dichotomique de l'energie du niveau fondamental********************************/
for(j=0;j<100;j++)
{
for(k=0;k<lg-1;k++)
{
v[k+1]=v[k]+h*(-((2*m)/(C*C)))*(((E+U)/2)-V((double)k))*x[k]; //Je test la convergence pour une valeur
x[k+1]=x[k]+h*v[k+1]; //comprise entre E et U
}
if (x[(int)((l/h))] >0) //Si cette convergence ce fait de maniere identique a celle de E trouvée
//precedement alors c'est une energie plus proche de la valeur exacte
E=(E+U)/2;
if (x[(int)((l/h))]<0) //En revenche si la convergence ce fait de manire differente alors
//je peux augmenter ma borne inferieure (U) car mon energie n'est
U=(E+U)/2; //pas dans le demi intervalle inferieur de [E;U]
}
}
if (i==p )
{
if (q%2==0)
{
p=p+2;
printf("Energie=%lf\n",E); //Cette partie du code sert a distinguer les moment ou l'on doit
} //rechercher une fonction paire de ceux ou l'on recherche une
else //fontion impaire (depend de la valeur de n).
{
p=p+1;
printf("Energie=%lf\n",E);
}
q++;
n=n+1;
}
}
/**************************************Dessin de fonction d'onde correspondante*************************************/
n=b;
if (n%2==0)
parite=1; //On impose la parite de la fonction d'onde suivant le mode que l'on recherche
else
parite=0;
if (parite == 1)
{
x[0]=1.0;
v[0]=0.0;
}
else
{
x[0]=0.0;
v[0]=1.0;
}
for(i=0;i<lg-1;i++)
{
fprintf(pf1,"%lf \t %lf \t %lf \n",h*i,x[i],V(i)); //On imprime les resultats dans un fichier et on complete suivant la paritee
if (parite == 1)
fprintf(pf1,"%lf \t %lf \t %lf \n",-h*i,x[i],V(i));
if (parite != 1)
fprintf(pf1,"%lf \t %lf \t %lf \n",-h*i,-x[i],V(i));
v[i+1]=v[i]+h*(-((2*m)/(C*C)))*(E-V((double)i))*x[i];
x[i+1]=x[i]+h*v[i+1];
if (x[i]< -5) //il ne sert a rien d'imprimer plus etant donne que une fois que la fonction a
break; //commancer a diverger elle ne va pas s'arreter
if (x[i] >5)
break;
}
printf("\n\n");
close(pf1);
}
Conclusion :
Je tiens a remercier specialement F.Hebert qui a coodonner ce projet.
J'essairai de donner une version plus complete de la facon dont marche ce programme ainsi que les ameliorations a y apporter ulterieurement.
Desole pour les maniaques de l'orthographe mais moi c'est pas mon truc alors....
Codes Sources
A voir également
- Equation de schrodinger : calcul des états propres dans un puits (méthode type euler)
- Equation de droite ✓ - Forum VB.NET
- Equation arc de cercle - Forum Visual Basic 6
- Equation cercle - Forum C++ & C++ .NET
- Equation ellipse ✓ - Forum Delphi / Pascal
- Solveur d'équation ✓ - Forum Visual Basic 6
Vous n'êtes pas encore membre ?
inscrivez-vous, c'est gratuit et ça prend moins d'une minute !
Les membres obtiennent plus de réponses que les utilisateurs anonymes.
Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.
Le fait d'être membre vous permet d'avoir des options supplémentaires.