ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/UseTheForce/DarkSide/suttonchen.F90
(Generate patch)

Comparing trunk/OOPSE-4/src/UseTheForce/DarkSide/suttonchen.F90 (file contents):
Revision 2530 by chuckv, Fri Dec 30 00:18:28 2005 UTC vs.
Revision 2541 by chuckv, Tue Jan 10 21:14:32 2006 UTC

# Line 1 | Line 1
1 < !!
1 > !!
2   !! Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   !!
4   !! The University of Notre Dame grants you ("Licensee") a
# Line 66 | Line 66 | module suttonchen
66    logical, save :: haveRcut = .false.
67    logical, save :: haveMixingMap = .false.
68    logical, save :: useGeometricDistanceMixing = .false.
69 <
70 <
69 >  logical, save :: cleanArrays = .true.
70 >  logical, save :: arraysAllocated = .false.
71  
72  
73    character(len = statusMsgSize) :: errMesg
# Line 76 | Line 76 | module suttonchen
76    character(len = 200) :: errMsg
77    character(len=*), parameter :: RoutineName =  "Sutton-Chen MODULE"
78    !! Logical that determines if eam arrays should be zeroed
79  logical :: cleanme = .true.
79    logical :: nmflag  = .false.
80  
81  
# Line 273 | Line 272 | contains
272            MixingMap(i,j)%rcut = 2.0_dp *MixingMap(i,j)%alpha
273            MixingMap(i,j)%vpair_pot = MixingMap(i,j)%epsilon* &
274                 (MixingMap(i,j)%alpha/MixingMap(i,j)%rcut)**MixingMap(i,j)%n
275 +
276            if (i.ne.j) then
277               MixingMap(j,i)%epsilon    = MixingMap(i,j)%epsilon
278               MixingMap(j,i)%m          = MixingMap(i,j)%m
# Line 291 | Line 291 | contains
291  
292    !! routine checks to see if array is allocated, deallocates array if allocated
293    !! and then creates the array to the required size
294 <  subroutine allocateSC(status)
295 <    integer, intent(out) :: status
294 >  subroutine allocateSC()
295 >    integer :: status
296  
297   #ifdef IS_MPI
298      integer :: nAtomsInRow
# Line 300 | Line 300 | contains
300   #endif
301      integer :: alloc_stat
302  
303 <
303 >    
304      status = 0
305   #ifdef IS_MPI
306      nAtomsInRow = getNatomsInRow(plan_atom_row)
# Line 313 | Line 313 | contains
313      allocate(frho(nlocal),stat=alloc_stat)
314      if (alloc_stat /= 0) then
315         status = -1
316       return
316      end if
317  
318      if (allocated(rho)) deallocate(rho)
319      allocate(rho(nlocal),stat=alloc_stat)
320      if (alloc_stat /= 0) then
321         status = -1
323       return
322      end if
323  
324      if (allocated(dfrhodrho)) deallocate(dfrhodrho)
325      allocate(dfrhodrho(nlocal),stat=alloc_stat)
326      if (alloc_stat /= 0) then
327         status = -1
330       return
328      end if
329  
330   #ifdef IS_MPI
# Line 336 | Line 333 | contains
333      allocate(rho_tmp(nlocal),stat=alloc_stat)
334      if (alloc_stat /= 0) then
335         status = -1
339       return
336      end if
337  
338  
# Line 344 | Line 340 | contains
340      allocate(frho_row(nAtomsInRow),stat=alloc_stat)
341      if (alloc_stat /= 0) then
342         status = -1
347       return
343      end if
344      if (allocated(rho_row)) deallocate(rho_row)
345      allocate(rho_row(nAtomsInRow),stat=alloc_stat)
346      if (alloc_stat /= 0) then
347         status = -1
353       return
348      end if
349      if (allocated(dfrhodrho_row)) deallocate(dfrhodrho_row)
350      allocate(dfrhodrho_row(nAtomsInRow),stat=alloc_stat)
351      if (alloc_stat /= 0) then
352         status = -1
359       return
353      end if
354  
355  
# Line 366 | Line 359 | contains
359      allocate(frho_col(nAtomsInCol),stat=alloc_stat)
360      if (alloc_stat /= 0) then
361         status = -1
369       return
362      end if
363      if (allocated(rho_col)) deallocate(rho_col)
364      allocate(rho_col(nAtomsInCol),stat=alloc_stat)
365      if (alloc_stat /= 0) then
366         status = -1
375       return
367      end if
368      if (allocated(dfrhodrho_col)) deallocate(dfrhodrho_col)
369      allocate(dfrhodrho_col(nAtomsInCol),stat=alloc_stat)
370      if (alloc_stat /= 0) then
371         status = -1
381       return
372      end if
373  
374   #endif
375 <
375 >    if (status == -1) then
376 >       call handleError("SuttonChen:allocateSC","Error in allocating SC arrays")
377 >    end if
378 >    arraysAllocated = .true.
379    end subroutine allocateSC
380  
381    !! C sets rcut to be the largest cutoff of any atype
# Line 397 | Line 390 | contains
390  
391    end subroutine setCutoffSC
392  
393 + !! This array allocates module arrays if needed and builds mixing map.
394    subroutine clean_SC()
395 <
395 >    if (.not.arraysAllocated) call allocateSC()
396 >    if (.not.haveMixingMap) call createMixingMap()
397      ! clean non-IS_MPI first
398      frho = 0.0_dp
399      rho  = 0.0_dp
# Line 439 | Line 434 | contains
434  
435      ! check to see if we need to be cleaned at the start of a force loop
436  
437 +    if (cleanArrays) call clean_SC()
438 +    cleanArrays = .false.
439  
443
444
440   #ifdef IS_MPI
441      Atid1 = Atid_row(Atom1)
442      Atid2 = Atid_col(Atom2)
# Line 476 | Line 471 | contains
471      real(kind=dp) :: pot
472      integer :: i,j
473      integer :: atom
479    real(kind=dp) :: U,U1,U2
474      integer :: atype1
475      integer :: atid1
476      integer :: myid
477  
478  
485    cleanme = .true.
479      !! Scatter the electron density from  pre-pair calculation back to local atoms
480   #ifdef IS_MPI
481      call scatter(rho_row,rho,plan_atom_row,sc_err)
# Line 506 | Line 499 | contains
499      do atom = 1, nlocal
500         Myid = SCList%atidtoSctype(Atid(atom))
501         frho(atom) = -SCList%SCTypes(Myid)%c * &
502 <            SCList%SCTypes(Myid)%epsilon * sqrt(rho(i))
502 >            SCList%SCTypes(Myid)%epsilon * sqrt(rho(atom))
503  
504         dfrhodrho(atom) = 0.5_dp*frho(atom)/rho(atom)
505 <       pot = pot + u
505 >
506 >       pot = pot + frho(atom)
507      enddo
508  
509   #ifdef IS_MPI
# Line 553 | Line 547 | contains
547      real( kind = dp ) :: dvpdr
548      real( kind = dp ) :: drhodr
549      real( kind = dp ) :: dudr
556    real( kind = dp ) :: rcij
550      real( kind = dp ) :: drhoidr,drhojdr
551      real( kind = dp ) :: Fx,Fy,Fz
552 <    real( kind = dp ) :: r,d2pha,phb,d2phb
552 >    real( kind = dp ) :: d2pha,phb,d2phb
553      real( kind = dp ) :: pot_temp,vptmp
554      real( kind = dp ) :: epsilonij,aij,nij,mij,vcij
555      integer :: id1,id2
# Line 566 | Line 559 | contains
559      !Local Variables
560  
561      ! write(*,*) "Frho: ", Frho(atom1)
562 +    
563 +    cleanArrays = .true.
564  
570
565      dvpdr = 0.0E0_DP
566  
567  
# Line 593 | Line 587 | contains
587         mij       = MixingMap(mytype_atom1,mytype_atom2)%m
588         vcij      = MixingMap(mytype_atom1,mytype_atom2)%vpair_pot              
589  
590 <       vptmp = epsilonij*((aij/r)**nij)
590 >       vptmp = epsilonij*((aij/rij)**nij)
591  
592  
593 <       dvpdr = -nij*vptmp/r
594 <       drhodr = -mij*((aij/r)**mij)/r
593 >       dvpdr = -nij*vptmp/rij
594 >       drhodr = -mij*((aij/rij)**mij)/rij
595  
596        
597         dudr = drhodr*(dfrhodrho(atom1)+dfrhodrho(atom2)) &

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines