Revision 10 src/egrad_siesta.f90

egrad_siesta.f90 (revision 10)
3 3
  ! This routines calculates the energy and the gradient of 
4 4
  ! a molecule, using Siesta 
5 5

  
6
  use Path_module, only : Nat, renum,Order,OrderInv, Coord, ProgExe
6
  use Path_module, only : Nat, renum,Order,OrderInv, Coord, ProgExe,v_direct,latr
7 7
  use Io_module
8 8

  
9 9
  IMPLICIT NONE
......
73 73
  logical           :: debug, TChk
74 74
  LOGICAL, SAVE :: first=.true.
75 75

  
76
  REAL(KREAL) :: Pi
76
  REAL(KREAL), ALLOCATABLE :: X(:),Y(:),Z(:) ! Nat
77
  REAL(KREAL), ALLOCATABLE :: GeomCartLoc(:,:) !Nat,3
77 78

  
78
  INTEGER(KINT) :: iat, i, n3at,Idx
79
  INTEGER(KINT) :: iat, i, J,Idx
79 80

  
80 81
  CHARACTER(LCHARS) :: FileIn, FileOut, FileForces,FileE
81 82

  
......
84 85

  
85 86
  ! ======================================================================
86 87

  
87
  Pi=dacos(-1.0d0)
88
  n3at=3*nat
89

  
90 88
  debug=valid('EGRAD').OR.valid('egrad_siesta')
91 89

  
92 90
  if (debug) Call Header("Entering Egrad_siesta")
......
103 101
     END IF
104 102
  END IF
105 103

  
104
  ALLOCATE(GeomCartLoc(Nat,3))
105

  
106 106
  gradcart=0.
107 107
  E=0.
108 108
  FileIn=Trim(CalcName) // "." // Trim(ISuffix)
......
122 122

  
123 123
  WRITE(IOTMP,*) 
124 124

  
125
  IF (V_Direct=='DIRECT') THEN
126
     ALLOCATE(X(Nat),Y(Nat),Z(Nat))
127
     X=0.
128
     Y=0.
129
     Z=0.
130
     DO I=1,Nat
131
        DO J=1,3
132
           X(I)=X(I)+GeomCart(I,J)*Latr(1,J)
133
           Y(I)=Y(I)+GeomCart(I,J)*Latr(2,J)
134
           Z(I)=Z(I)+GeomCart(I,J)*Latr(3,J)
135
        END DO
136
     END DO
137
     GeomCartLoc(1:Nat,1)=X(1:Nat)
138
     GeomCartLoc(1:Nat,2)=Y(1:Nat)
139
     GeomCartLoc(1:Nat,3)=Z(1:Nat)
140
     DEALLOCATE(X,Y,Z)
141
  ELSE
142
     GeomCartLoc=GeomCart
143
  END IF
144

  
125 145
  WRITE(IOTMP,'(1X,A)')  '%block AtomicCoordinatesAndAtomicSpecies' 
126 146

  
127 147
  DO I=1,Nat
......
206 226
     WRITE(*,'(1X,3(1X,F15.8))') GradCart(3*I-2:3*I)
207 227
  END DO
208 228

  
229
  DeAllocate(GeomCartLoc)
230

  
209 231
!  Call Die('Egrad_Siesta','Ah ah.',UNIT=IOOUT)
210 232
  if (debug) Call header("Egrad_siesta Over")
211 233

  
234

  
235

  
212 236
  RETURN
213 237

  
214 238
  ! ======================================================================

Also available in: Unified diff