Discussion :

[romanticide]

Bonjour,

Est-ce que qq un sait comment programmer en fortran l'énergie cinétique d'un système de particule ?

Merci !

[Ladgalen]

Tu veux juste la calculer ?

Il suffit que tu connaisses la vitesse de chacun de tes atomes ou de tes particules, ça n'a rien de mystérieux :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
 
integer, parameter                 :: natom = 1000
integer                            :: iat
double precision, dimension(nat,3) :: v
double precision, dimension(nat)   :: masse
double precision                   :: ecinet
 
ecinet = 0.d0
do iat = 1, nat
  ecinet = ecinet + 0.5.d0 * masse(iat) * (v(iat,1)**2 + v(iat,2)**2 + v(iat,3)**2)
  ! ou alors ecinet = ecinet + 0.5 * masse(iat) * sum( v(iat,:)**2 )
end do

Ce n'est qu'une simple boucle !

[romanticide]

je connais la vitesse selon la direction x, selon la direction y et selon la direction z
Est-ce que je dois faire le calcul de l'énergie pour chaque direction ?
je ne comprends pas comment obtenir l'énergie du système si tout est décomposé selon les 3 dimensions

merci pour le code ça va m'aider

[Ladgalen]

Je pense que tu as plus un problème de math que de fortran. Tout est une affaire de vecteur, donc la seule chose que tu connaissent ce sont les coordonnées de ces vecteurs. Vecteur vitesse, vecteur position de tes particules ... Ensuite à toi de calculer la norme du vecteur pour calculer la vitesse, ou la distance ect ...

L'énergie cinétique c'est E = 1/2 m v^2

Tu sais que la vitesse c'est un vecteur, donc tu connais logiquement les coordonnées du vecteur vitesse : vx, vy, vz. Ce que tu appelles la vitesse selon la direction x y ou z. Tu peux les représenter en fortran par un tableau de dimension 3

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

double precision, dimension(3) :: v

Pour obtenir la vitesse de ta particule, disons la valeur de la vitesse, elle correspond à la norme du vecteur vitesse donc :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

||v|| = sqrt( vx^2 + vy^2 + vz^2 )

C'est ce nombre qui intervient dans ton énergie cinétique. Dans mon code j'ai défini un tableau vitesse qui contient toutes les vitesses de toutes les particules

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

double precision, dimension(nat,3) :: v

donc

v(392,1) = vx de la particule 382
v(392,2) = vy de la particule 382
v(392,3) = vz de la particule 382

De la tu peux en déduire l'énergie cinétique de ta particule 392.

Ecinet de 392 = 0.5 * masse de 392 * (v(392,1)**2 + v(392,2)**2 + v(392,3)**2)

Est-ce que je dois faire le calcul de l'énergie pour chaque direction ?

ça n'a pas d'importance, tu peux définir une énergie cinétique selon x selon y et selon z et faire la somme des trois pour avoir l'énergie cinétique ou calculer directement la somme

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
integer, parameter                 :: natom = 1000
integer                            :: iat
double precision, dimension(nat,3) :: v
double precision, dimension(nat)   :: masse
double precision                   :: ecinet, ecinetx, ecinety, ecinetz
double precision                   :: eci
double precision                   :: vx, vy, vz
 
! energie cinetique de la particule 392
vx = v(392,1)
vy = v(392,2)
vz = v(392,3)
eci = 0.5 * masse(392) * ( vx**2 + vy**2 + vz**2 )
 
! que tu peux ecrire 
eci = 0.5 * masse(392) * vx**2 + 0.5 * masse(392) * vy**2 + 0.5 * masse(392) * vz**2 
 
! energie cinetique totale en calculant selon chaque direction
ecinet = 0.d0
do iat = 1, nat
  vx = v(iat,1)
  vy = v(iat,2)
  vz = v(iat,3)
 
  ecinetx = 0.5 * masse(iat) * vx**2
  ecinety = 0.5 * masse(iat) * vy**2
  ecinetz = 0.5 * masse(iat) * vz**2
 
  ecinet = ecinet + ecinetx + ecinety + ecinetz
end do
 
! energie cinetique totale d'un coup (c'est juste une factorisation de 0.5*masse)
ecinet = 0.d0
do iat = 1, nat
  vx = v(iat,1)
  vy = v(iat,2)
  vz = v(iat,3)
  ecinet = ecinet + 0.5.d0 * masse(iat) * (vx**2 + vy**2 + vz**2)
end do

Si ça tu veux, sur ma page web tu peux récupérer un code complet de dynamique moléculaire. Il est assez petit mais il contient à peu près toutes les bases. Pour des raisons diverses et variés il est écrit en C ... Tu peux notamment regarder le fichier moyenne.c dans lequel je calcul l'énergie cinétique la température la pression ect ...

[romanticide]

merci BEAUCOUP !!!