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

Comparing trunk/OOPSE-4/src/UseTheForce/DarkSide/LJ.F90 (file contents):
Revision 2361 by gezelter, Wed Oct 12 21:00:50 2005 UTC vs.
Revision 2461 by gezelter, Mon Nov 21 22:59:02 2005 UTC

# Line 43 | Line 43
43   !! Calculates Long Range forces Lennard-Jones interactions.
44   !! @author Charles F. Vardeman II
45   !! @author Matthew Meineke
46 < !! @version $Id: LJ.F90,v 1.18 2005-10-12 21:00:49 gezelter Exp $, $Date: 2005-10-12 21:00:49 $, $Name: not supported by cvs2svn $, $Revision: 1.18 $
46 > !! @version $Id: LJ.F90,v 1.19 2005-11-21 22:59:01 gezelter Exp $, $Date: 2005-11-21 22:59:01 $, $Name: not supported by cvs2svn $, $Revision: 1.19 $
47  
48  
49   module lj
# Line 102 | Line 102 | module lj
102  
103    public :: newLJtype
104    public :: setLJDefaultCutoff
105  public :: setLJUniformCutoff
106  public :: setLJCutoffByTypes
105    public :: getSigma
106    public :: getEpsilon
107    public :: useGeometricMixing
# Line 164 | Line 162 | contains
162         call createMixingMap()
163      end if
164    end subroutine setLJDefaultCutoff
167
168  subroutine setLJUniformCutoff(thisRcut, shiftedPot)
169    real(kind=dp), intent(in) :: thisRcut
170    logical, intent(in) :: shiftedPot
171    integer :: nLJtypes, i, j
172
173    if (LJMap%currentLJtype == 0) then
174       call handleError("LJ", "No members in LJMap")
175       return
176    end if
177
178    nLJtypes = LJMap%nLJtypes
179    if (.not. allocated(MixingMap)) then
180       allocate(MixingMap(nLJtypes, nLJtypes))
181    endif
182
183    do i = 1, nLJtypes
184       do j = 1, nLJtypes
185          MixingMap(i,j)%rCut = thisRcut
186          MixingMap(i,j)%shiftedPot = shiftedPot
187          MixingMap(i,j)%rCutWasSet = .true.
188       enddo
189    enddo
190    call createMixingMap()
191  end subroutine setLJUniformCutoff
192
193  subroutine setLJCutoffByTypes(atid1, atid2, thisRcut, shiftedPot)
194    integer, intent(in) :: atid1, atid2
195    real(kind=dp), intent(in) :: thisRcut
196    logical, intent(in) :: shiftedPot
197    integer :: nLJtypes, ljt1, ljt2
198
199    if (LJMap%currentLJtype == 0) then
200       call handleError("LJ", "No members in LJMap")
201       return
202    end if
203
204    nLJtypes = LJMap%nLJtypes
205    if (.not. allocated(MixingMap)) then
206       allocate(MixingMap(nLJtypes, nLJtypes))
207    endif
208
209    ljt1 = LJMap%atidToLJtype(atid1)
210    ljt2 = LJMap%atidToLJtype(atid2)
211
212    MixingMap(ljt1,ljt2)%rCut = thisRcut
213    MixingMap(ljt1,ljt2)%shiftedPot = shiftedPot
214    MixingMap(ljt1,ljt2)%rCutWasSet = .true.
215    MixingMap(ljt2,ljt1)%rCut = thisRcut
216    MixingMap(ljt2,ljt1)%shiftedPot = shiftedPot
217    MixingMap(ljt2,ljt1)%rCutWasSet = .true.
165  
219    call createMixingMap()
220  end subroutine setLJCutoffByTypes
221
166    function getSigma(atid) result (s)
167      integer, intent(in) :: atid
168      integer :: ljt1
# Line 297 | Line 241 | contains
241  
242            MixingMap(i,j)%sigma6 = (MixingMap(i,j)%sigma)**6
243  
244 <          if (MixingMap(i,j)%rCutWasSet) then
245 <             rcut6 = (MixingMap(i,j)%rcut)**6
244 >          if (haveDefaultCutoff) then
245 >             rcut6 = defaultCutoff**6
246 >             tp6    = MixingMap(i,j)%sigma6/rcut6
247 >             tp12    = tp6**2          
248 >             MixingMap(i,j)%delta =-4.0_DP*MixingMap(i,j)%epsilon*(tp12 - tp6)
249 >             MixingMap(i,j)%shiftedPot = defaultShift
250            else
251 <             if (haveDefaultCutoff) then
252 <                rcut6 = defaultCutoff**6
253 <                doShift = defaultShift
306 <             else
307 <                call handleError("LJ", "No specified or default cutoff value!")
308 <             endif
309 <          endif
310 <          
311 <          tp6    = MixingMap(i,j)%sigma6/rcut6
312 <          tp12    = tp6**2          
313 <          MixingMap(i,j)%delta =-4.0_DP*MixingMap(i,j)%epsilon*(tp12 - tp6)
314 <          MixingMap(i,j)%shiftedPot = doShift
251 >             MixingMap(i,j)%delta = 0.0_DP
252 >             MixingMap(i,j)%shiftedPot = defaultShift
253 >          endif          
254  
255            if (i.ne.j) then
256               MixingMap(j,i)%sigma      = MixingMap(i,j)%sigma
# Line 331 | Line 270 | contains
270      
271    end subroutine createMixingMap
272    
273 <  subroutine do_lj_pair(atom1, atom2, d, rij, r2, sw, vpair, fpair, &
273 >  subroutine do_lj_pair(atom1, atom2, d, rij, r2, rcut, sw, vpair, fpair, &
274         pot, f, do_pot)
275 <
275 >    
276      integer, intent(in) ::  atom1, atom2
277      integer :: atid1, atid2, ljt1, ljt2
278 <    real( kind = dp ), intent(in) :: rij, r2
278 >    real( kind = dp ), intent(in) :: rij, r2, rcut
279      real( kind = dp ) :: pot, sw, vpair
280      real( kind = dp ), dimension(3,nLocal) :: f    
281      real( kind = dp ), intent(in), dimension(3) :: d
# Line 349 | Line 288 | contains
288      real( kind = dp ) :: pot_temp, dudr
289      real( kind = dp ) :: sigma6
290      real( kind = dp ) :: epsilon
291 <    real( kind = dp ) :: r6
292 <    real( kind = dp ) :: t6
293 <    real( kind = dp ) :: t12
291 >    real( kind = dp ) :: r6, rc6
292 >    real( kind = dp ) :: t6, tp6
293 >    real( kind = dp ) :: t12, tp12
294      real( kind = dp ) :: delta
295      logical :: isSoftCore, shiftedPot
296      integer :: id1, id2, localError
# Line 387 | Line 326 | contains
326        
327         pot_temp = 4.0E0_DP * epsilon * t6
328         if (shiftedPot) then
329 +          rc6 = rcut**6
330 +          tp6 = sigma6 / rc6
331 +          delta =-4.0_DP*epsilon*(tp6)
332            pot_temp = pot_temp + delta
333         endif
334        
# Line 397 | Line 339 | contains
339      else
340         pot_temp = 4.0E0_DP * epsilon * (t12 - t6)
341         if (shiftedPot) then
342 +          rc6 = rcut**6
343 +          tp6 = sigma6 / rc6
344 +          tp12 = tp6*tp6
345 +          delta =-4.0_DP*epsilon*(tp12 - tp6)
346            pot_temp = pot_temp + delta
347         endif
348        

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines