1c
2c
3c     ###################################################
4c     ##  COPYRIGHT (C)  1993  by  Jay William Ponder  ##
5c     ##              All Rights Reserved              ##
6c     ###################################################
7c
8c     ###########################################################
9c     ##                                                       ##
10c     ##  subroutine eurey  --  Urey-Bradley potential energy  ##
11c     ##                                                       ##
12c     ###########################################################
13c
14c
15c     "eurey" calculates the Urey-Bradley 1-3 interaction energy
16c
17c
18      subroutine eurey
19      use atoms
20      use bound
21      use energi
22      use group
23      use urey
24      use urypot
25      use usage
26      implicit none
27      integer i,ia,ic
28      real*8 e,ideal,force
29      real*8 dt,dt2,fgrp
30      real*8 xac,yac,zac,rac
31      logical proceed
32c
33c
34c     zero out the Urey-Bradley interaction energy
35c
36      eub = 0.0d0
37      if (nurey .eq. 0)  return
38c
39c     OpenMP directives for the major loop structure
40c
41!$OMP PARALLEL default(private) shared(nurey,iury,ul,uk,
42!$OMP& use,x,y,z,cury,qury,ureyunit,use_group,use_polymer)
43!$OMP& shared(eub)
44!$OMP DO reduction(+:eub) schedule(guided)
45c
46c     calculate the Urey-Bradley 1-3 energy term
47c
48      do i = 1, nurey
49         ia = iury(1,i)
50         ic = iury(3,i)
51         ideal = ul(i)
52         force = uk(i)
53c
54c     decide whether to compute the current interaction
55c
56         proceed = .true.
57         if (use_group)  call groups (proceed,fgrp,ia,ic,0,0,0,0)
58         if (proceed)  proceed = (use(ia) .or. use(ic))
59c
60c     compute the value of the 1-3 distance deviation
61c
62         if (proceed) then
63            xac = x(ia) - x(ic)
64            yac = y(ia) - y(ic)
65            zac = z(ia) - z(ic)
66            if (use_polymer)  call image (xac,yac,zac)
67            rac = sqrt(xac*xac + yac*yac + zac*zac)
68            dt = rac - ideal
69            dt2 = dt * dt
70c
71c     calculate the Urey-Bradley energy for this interaction
72c
73            e = ureyunit * force * dt2 * (1.0d0+cury*dt+qury*dt2)
74c
75c     scale the interaction based on its group membership
76c
77            if (use_group)  e = e * fgrp
78c
79c     increment the total Urey-Bradley energy
80c
81            eub = eub + e
82         end if
83      end do
84c
85c     OpenMP directives for the major loop structure
86c
87!$OMP END DO
88!$OMP END PARALLEL
89      return
90      end
91