1c 2c 3c ################################################### 4c ## COPYRIGHT (C) 1993 by Jay William Ponder ## 5c ## All Rights Reserved ## 6c ################################################### 7c 8c ############################################################# 9c ## ## 10c ## subroutine eurey3 -- Urey-Bradley energy & analysis ## 11c ## ## 12c ############################################################# 13c 14c 15c "eurey3" calculates the Urey-Bradley energy; also 16c partitions the energy among the atoms 17c 18c 19 subroutine eurey3 20 use action 21 use analyz 22 use atomid 23 use atoms 24 use bound 25 use energi 26 use group 27 use inform 28 use iounit 29 use urey 30 use urypot 31 use usage 32 implicit none 33 integer i,ia,ib,ic 34 real*8 e,ideal,force 35 real*8 dt,dt2,fgrp 36 real*8 xac,yac,zac,rac 37 logical proceed 38 logical header,huge 39c 40c 41c zero out the Urey-Bradley energy and partitioning terms 42c 43 neub = 0 44 eub = 0.0d0 45 do i = 1, n 46 aeub(i) = 0.0d0 47 end do 48 if (nurey .eq. 0) return 49c 50c print header information if debug output was requested 51c 52 header = .true. 53 if (debug .and. nurey.ne.0) then 54 header = .false. 55 write (iout,10) 56 10 format (/,' Individual Urey-Bradley Interactions :', 57 & //,' Type',18x,'Atom Names',18x,'Ideal', 58 & 4x,'Actual',6x,'Energy',/) 59 end if 60c 61c OpenMP directives for the major loop structure 62c 63!$OMP PARALLEL default(private) shared(nurey,iury,ul,uk, 64!$OMP& use,x,y,z,cury,qury,ureyunit,use_group,use_polymer, 65!$OMP& name,verbose,debug,header,iout) 66!$OMP& shared(eub,neub,aeub) 67!$OMP DO reduction(+:eub,neub,aeub) schedule(guided) 68c 69c calculate the Urey-Bradley 1-3 energy term 70c 71 do i = 1, nurey 72 ia = iury(1,i) 73 ib = iury(2,i) 74 ic = iury(3,i) 75 ideal = ul(i) 76 force = uk(i) 77c 78c decide whether to compute the current interaction 79c 80 proceed = .true. 81 if (use_group) call groups (proceed,fgrp,ia,ic,0,0,0,0) 82 if (proceed) proceed = (use(ia) .or. use(ic)) 83c 84c compute the value of the 1-3 distance deviation 85c 86 if (proceed) then 87 xac = x(ia) - x(ic) 88 yac = y(ia) - y(ic) 89 zac = z(ia) - z(ic) 90 if (use_polymer) call image (xac,yac,zac) 91 rac = sqrt(xac*xac + yac*yac + zac*zac) 92 dt = rac - ideal 93 dt2 = dt * dt 94c 95c calculate the Urey-Bradley energy for this interaction 96c 97 e = ureyunit * force * dt2 * (1.0d0+cury*dt+qury*dt2) 98c 99c scale the interaction based on its group membership 100c 101 if (use_group) e = e * fgrp 102c 103c increment the total Urey-Bradley energy 104c 105 neub = neub + 1 106 eub = eub + e 107 aeub(ia) = aeub(ia) + 0.5d0*e 108 aeub(ic) = aeub(ic) + 0.5d0*e 109c 110c print a message if the energy of this interaction is large 111c 112 huge = (e .gt. 5.0d0) 113 if (debug .or. (verbose.and.huge)) then 114 if (header) then 115 header = .false. 116 write (iout,20) 117 20 format (/,' Individual Urey-Bradley Interactions :', 118 & //,' Type',18x,'Atom Names',18x,'Ideal', 119 & 4x,'Actual',6x,'Energy',/) 120 end if 121 write (iout,30) ia,name(ia),ib,name(ib), 122 & ic,name(ic),ideal,rac,e 123 30 format (' UreyBrad',2x,i7,'-',a3,i7,'-',a3, 124 & i7,'-',a3,2x,2f10.4,f12.4) 125 end if 126 end if 127 end do 128c 129c OpenMP directives for the major loop structure 130c 131!$OMP END DO 132!$OMP END PARALLEL 133 return 134 end 135