root / src / ConvertZmat_cart_3.f90 @ 2
Historique | Voir | Annoter | Télécharger (2,43 ko)
| 1 |
!C================================================================ |
|---|---|
| 2 |
!C Converti les positions Zmat en coordonnes cartesiennes |
| 3 |
!C This version deals with the third atoms in A B C |
| 4 |
!C Used in Mix to Cart. |
| 5 |
!C================================================================ |
| 6 |
|
| 7 |
SUBROUTINE ConvertZmat_cart_3(iat,ind_zmat,d,a_val,x,y,z) |
| 8 |
|
| 9 |
use Path_module, only : Nat, KINT, KREAL |
| 10 |
|
| 11 |
IMPLICIT NONE |
| 12 |
|
| 13 |
integer(KINT) :: iat,n1,n2,n3 |
| 14 |
real(KREAL) :: x(nat),y(nat),z(nat) |
| 15 |
real(KREAL) :: d,a_val |
| 16 |
integer(KINT) :: ind_zmat(Nat,5) |
| 17 |
|
| 18 |
real(KREAL) :: vx1,vy1,vz1,norm1 |
| 19 |
real(KREAL) :: vvx1,vvy1,vvz1,normv1 |
| 20 |
real(KREAL) :: vx2,vy2,vz2,norm2 |
| 21 |
real(KREAL) :: vvx2,vvy2,vvz2,normv2 |
| 22 |
real(KREAL) :: vx3,vy3,vz3,norm3 |
| 23 |
real(KREAL) :: vx4,vy4,vz4,norm4 |
| 24 |
real(KREAL) :: a11_z1,a12_z1 |
| 25 |
real(KREAL) :: a11_z2,a12_z2 |
| 26 |
real(KREAL) :: a11_y,a12_y |
| 27 |
|
| 28 |
|
| 29 |
|
| 30 |
!C ind_zmat(1) contient le numero de l'atome, (2) celui par rapport auquel on definit la distance... |
| 31 |
|
| 32 |
n1=ind_zmat(iat,2) |
| 33 |
n2=ind_zmat(iat,3) |
| 34 |
|
| 35 |
! WRITE(*,*) iat,n1,n2 |
| 36 |
|
| 37 |
! WRITE(*,*) "d,val,di",d,a_val,a_dih,a_val*180./3.141592654 |
| 38 |
|
| 39 |
CALL vecteur(n1,n2,x,y,z,vx1,vy1,vz1,norm1) |
| 40 |
|
| 41 |
vvx1=vx1 |
| 42 |
vvy1=vy1 |
| 43 |
vvz1=vz1 |
| 44 |
|
| 45 |
!c rotation autour de z de v1 de phi (a11_z1=cos(phi) et a12_z1=sin(phi) ) tq |
| 46 |
!C v1 soit dans le plan xOz |
| 47 |
|
| 48 |
normv1=dsqrt(vvx1*vvx1+vvy1*vvy1) |
| 49 |
|
| 50 |
IF (normv1 .GE. 1.D-6) THEN |
| 51 |
a11_z1 = vvx1/normv1 |
| 52 |
a12_z1 = vvy1/normv1 |
| 53 |
ELSE |
| 54 |
a11_z1 = 1 |
| 55 |
a12_z1 = 0 |
| 56 |
END IF |
| 57 |
|
| 58 |
CALL rota_z(vx1,vy1,vz1,a11_z1,-a12_z1) |
| 59 |
|
| 60 |
!c rotation autour de y de v1 de theta (a11_y=cos(theta) et a12_y=sin(theta) ) tq |
| 61 |
!C v1 soit sur l'axe Ox |
| 62 |
|
| 63 |
|
| 64 |
IF (norm1 .GE. 1.D-8) THEN |
| 65 |
a11_y = vz1/norm1 |
| 66 |
a12_y = vx1/norm1 |
| 67 |
ELSE |
| 68 |
a11_y = 1 |
| 69 |
a12_y = 0 |
| 70 |
END IF |
| 71 |
|
| 72 |
CALL rota_y(vx1,vy1,vz1,a11_y,a12_y) |
| 73 |
|
| 74 |
|
| 75 |
!c calcul le vecteur de l atome dans la nouvelle orientation |
| 76 |
|
| 77 |
vx4=d*sin(a_val) |
| 78 |
vy4=0. |
| 79 |
vz4=d*cos(a_val) |
| 80 |
|
| 81 |
!c calcul le vecteur de l atome |
| 82 |
!c on tourne en sens invers |
| 83 |
CALL rota_y(vx4,vy4,vz4,a11_y ,-a12_y ) |
| 84 |
CALL rota_z(vx4,vy4,vz4,a11_z1, a12_z1) |
| 85 |
|
| 86 |
!c calcul l atome a partir de v4 |
| 87 |
x(iat)=vx4+x(n1) |
| 88 |
y(iat)=vy4+y(n1) |
| 89 |
z(iat)=vz4+z(n1) |
| 90 |
|
| 91 |
END |
| 92 |
|