[ViewVC] Diff of: group/trunk/OOPSE/libmdtools/do

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents):
Revision 843 by mmeineke, Wed Oct 29 20:41:39 2003 UTC vs.
Revision 900 by chuckv, Tue Jan 6 18:54:57 2004 UTC

#	Line 4 \| Line 4
4
5		!! @author Charles F. Vardeman II
6		!! @author Matthew Meineke
7	<	!! @version $Id: do_Forces.F90,v 1.35 2003-10-29 20:41:39 mmeineke Exp $, $Date: 2003-10-29 20:41:39 $, $Name: not supported by cvs2svn $, $Revision: 1.35 $
7	>	!! @version $Id: do_Forces.F90,v 1.45 2004-01-06 18:54:57 chuckv Exp $, $Date: 2004-01-06 18:54:57 $, $Name: not supported by cvs2svn $, $Revision: 1.45 $
8
9		module do_Forces
10		use force_globals
#	Line 30 \| Line 30 \| module do_Forces
30		#define __FORTRAN90
31		#include "fForceField.h"
32
33	<	logical, save :: do_forces_initialized = .false., haveRlist = .false.
33	>	logical, save :: haveRlist = .false.
34	>	logical, save :: haveNeighborList = .false.
35		logical, save :: havePolicies = .false.
36	+	logical, save :: haveSIMvariables = .false.
37	+	logical, save :: havePropertyMap = .false.
38	+	logical, save :: haveSaneForceField = .false.
39		logical, save :: FF_uses_LJ
40		logical, save :: FF_uses_sticky
41		logical, save :: FF_uses_dipoles
42		logical, save :: FF_uses_RF
43		logical, save :: FF_uses_GB
44		logical, save :: FF_uses_EAM
45	+	logical, save :: SIM_uses_LJ
46	+	logical, save :: SIM_uses_sticky
47	+	logical, save :: SIM_uses_dipoles
48	+	logical, save :: SIM_uses_RF
49	+	logical, save :: SIM_uses_GB
50	+	logical, save :: SIM_uses_EAM
51	+	logical, save :: SIM_requires_postpair_calc
52	+	logical, save :: SIM_requires_prepair_calc
53	+	logical, save :: SIM_uses_directional_atoms
54	+	logical, save :: SIM_uses_PBC
55
56		real(kind=dp), save :: rlist, rlistsq
57
#	Line 45 \| Line 59 \| module do_Forces
59		public :: do_force_loop
60		public :: setRlistDF
61
62	+
63		#ifdef PROFILE
64	<	real(kind = dp) :: forceTime
65	<	real(kind = dp) :: forceTimeInitial, forceTimeFinal
66	<	real(kind = dp) :: globalForceTime
67	<	real(kind = dp) :: maxForceTime
53	<	integer, save :: nloops = 0
64	>	public :: getforcetime
65	>	real, save :: forceTime = 0
66	>	real :: forceTimeInitial, forceTimeFinal
67	>	integer :: nLoops
68		#endif
69
70	+	type :: Properties
71	+	logical :: is_lj = .false.
72	+	logical :: is_sticky = .false.
73	+	logical :: is_dp = .false.
74	+	logical :: is_gb = .false.
75	+	logical :: is_eam = .false.
76	+	real(kind=DP) :: dipole_moment = 0.0_DP
77	+	end type Properties
78	+
79	+	type(Properties), dimension(:),allocatable :: PropertyMap
80	+
81		contains
82
83		subroutine setRlistDF( this_rlist )
#	Line 63 \| Line 88 \| contains
88		rlistsq = rlist * rlist
89
90		haveRlist = .true.
66	–	if( havePolicies ) do_forces_initialized = .true.
91
92		end subroutine setRlistDF
93
94	+	subroutine createPropertyMap(status)
95	+	integer :: nAtypes
96	+	integer :: status
97	+	integer :: i
98	+	logical :: thisProperty
99	+	real (kind=DP) :: thisDPproperty
100	+
101	+	status = 0
102	+
103	+	nAtypes = getSize(atypes)
104	+
105	+	if (nAtypes == 0) then
106	+	status = -1
107	+	return
108	+	end if
109	+
110	+	if (.not. allocated(PropertyMap)) then
111	+	allocate(PropertyMap(nAtypes))
112	+	endif
113	+
114	+	do i = 1, nAtypes
115	+	call getElementProperty(atypes, i, "is_LJ", thisProperty)
116	+	PropertyMap(i)%is_LJ = thisProperty
117	+	call getElementProperty(atypes, i, "is_DP", thisProperty)
118	+	PropertyMap(i)%is_DP = thisProperty
119	+
120	+	if (thisProperty) then
121	+	call getElementProperty(atypes, i, "dipole_moment", thisDPproperty)
122	+	PropertyMap(i)%dipole_moment = thisDPproperty
123	+	endif
124	+
125	+	call getElementProperty(atypes, i, "is_Sticky", thisProperty)
126	+	PropertyMap(i)%is_Sticky = thisProperty
127	+	call getElementProperty(atypes, i, "is_GB", thisProperty)
128	+	PropertyMap(i)%is_GB = thisProperty
129	+	call getElementProperty(atypes, i, "is_EAM", thisProperty)
130	+	PropertyMap(i)%is_EAM = thisProperty
131	+	end do
132	+
133	+	havePropertyMap = .true.
134	+
135	+	end subroutine createPropertyMap
136	+
137	+	subroutine setSimVariables()
138	+	SIM_uses_LJ = SimUsesLJ()
139	+	SIM_uses_sticky = SimUsesSticky()
140	+	SIM_uses_dipoles = SimUsesDipoles()
141	+	SIM_uses_RF = SimUsesRF()
142	+	SIM_uses_GB = SimUsesGB()
143	+	SIM_uses_EAM = SimUsesEAM()
144	+	SIM_requires_postpair_calc = SimRequiresPostpairCalc()
145	+	SIM_requires_prepair_calc = SimRequiresPrepairCalc()
146	+	SIM_uses_directional_atoms = SimUsesDirectionalAtoms()
147	+	SIM_uses_PBC = SimUsesPBC()
148	+
149	+	haveSIMvariables = .true.
150	+
151	+	return
152	+	end subroutine setSimVariables
153	+
154	+	subroutine doReadyCheck(error)
155	+	integer, intent(out) :: error
156	+
157	+	integer :: myStatus
158	+
159	+	error = 0
160	+
161	+	if (.not. havePropertyMap) then
162	+
163	+	myStatus = 0
164	+
165	+	call createPropertyMap(myStatus)
166	+
167	+	if (myStatus .ne. 0) then
168	+	write(default_error, *) 'createPropertyMap failed in do_Forces!'
169	+	error = -1
170	+	return
171	+	endif
172	+	endif
173	+
174	+	if (.not. haveSIMvariables) then
175	+	call setSimVariables()
176	+	endif
177	+
178	+	if (.not. haveRlist) then
179	+	write(default_error, *) 'rList has not been set in do_Forces!'
180	+	error = -1
181	+	return
182	+	endif
183	+
184	+	if (SIM_uses_LJ .and. FF_uses_LJ) then
185	+	if (.not. havePolicies) then
186	+	write(default_error, *) 'LJ mixing Policies have not been set in do_Forces!'
187	+	error = -1
188	+	return
189	+	endif
190	+	endif
191	+
192	+	if (.not. haveNeighborList) then
193	+	write(default_error, *) 'neighbor list has not been initialized in do_Forces!'
194	+	error = -1
195	+	return
196	+	end if
197	+
198	+	if (.not. haveSaneForceField) then
199	+	write(default_error, *) 'Force Field is not sane in do_Forces!'
200	+	error = -1
201	+	return
202	+	end if
203	+
204	+	#ifdef IS_MPI
205	+	if (.not. isMPISimSet()) then
206	+	write(default_error,*) "ERROR: mpiSimulation has not been initialized!"
207	+	error = -1
208	+	return
209	+	endif
210	+	#endif
211	+	return
212	+	end subroutine doReadyCheck
213	+
214	+
215		subroutine init_FF(LJMIXPOLICY, use_RF_c, thisStat)
216
217		integer, intent(in) :: LJMIXPOLICY
#	Line 111 \| Line 256 \| contains
256		call getMatchingElementList(atypes, "is_EAM", .true., nMatches, MatchList)
257		if (nMatches .gt. 0) FF_uses_EAM = .true.
258
259	+	!! Assume sanity (for the sake of argument)
260	+	haveSaneForceField = .true.
261	+
262		!! check to make sure the FF_uses_RF setting makes sense
263
264		if (FF_uses_dipoles) then
#	Line 122 \| Line 270 \| contains
270		if (FF_uses_RF) then
271		write(default_error,*) 'Using Reaction Field with no dipoles? Huh?'
272		thisStat = -1
273	+	haveSaneForceField = .false.
274		return
275		endif
276		endif
#	Line 136 \| Line 285 \| contains
285		case default
286		write(default_error,*) 'unknown LJ Mixing Policy!'
287		thisStat = -1
288	+	haveSaneForceField = .false.
289		return
290		end select
291		if (my_status /= 0) then
292		thisStat = -1
293	+	haveSaneForceField = .false.
294		return
295		end if
296	+	havePolicies = .true.
297		endif
298
299		if (FF_uses_sticky) then
300		call check_sticky_FF(my_status)
301		if (my_status /= 0) then
302		thisStat = -1
303	+	haveSaneForceField = .false.
304		return
305		end if
306		endif
#	Line 156 \| Line 309 \| contains
309		if (FF_uses_EAM) then
310		call init_EAM_FF(my_status)
311		if (my_status /= 0) then
312	<	write(,) "init_EAM_FF returned a bad status"
312	>	write(default_error, *) "init_EAM_FF returned a bad status"
313		thisStat = -1
314	+	haveSaneForceField = .false.
315		return
316		end if
317		endif
318
165	–
166	–
319		if (FF_uses_GB) then
320		call check_gb_pair_FF(my_status)
321		if (my_status .ne. 0) then
322		thisStat = -1
323	+	haveSaneForceField = .false.
324		return
325		endif
326		endif
327
328		if (FF_uses_GB .and. FF_uses_LJ) then
329		endif
330	<	if (.not. do_forces_initialized) then
330	>	if (.not. haveNeighborList) then
331		!! Create neighbor lists
332	<	call expandNeighborList(getNlocal(), my_status)
332	>	call expandNeighborList(nLocal, my_status)
333		if (my_Status /= 0) then
334		write(default_error,*) "SimSetup: ExpandNeighborList returned error."
335		thisStat = -1
336		return
337		endif
338	+	haveNeighborList = .true.
339		endif
340
187	–
188	–	havePolicies = .true.
189	–	if( haveRlist ) do_forces_initialized = .true.
190	–
341		end subroutine init_FF
342
343
#	Line 196 \| Line 346 \| contains
346		subroutine do_force_loop(q, A, u_l, f, t, tau, pot, do_pot_c, do_stress_c, &
347		error)
348		!! Position array provided by C, dimensioned by getNlocal
349	<	real ( kind = dp ), dimension(3,getNlocal()) :: q
349	>	real ( kind = dp ), dimension(3,nLocal) :: q
350		!! Rotation Matrix for each long range particle in simulation.
351	<	real( kind = dp), dimension(9,getNlocal()) :: A
351	>	real( kind = dp), dimension(9,nLocal) :: A
352		!! Unit vectors for dipoles (lab frame)
353	<	real( kind = dp ), dimension(3,getNlocal()) :: u_l
353	>	real( kind = dp ), dimension(3,nLocal) :: u_l
354		!! Force array provided by C, dimensioned by getNlocal
355	<	real ( kind = dp ), dimension(3,getNlocal()) :: f
355	>	real ( kind = dp ), dimension(3,nLocal) :: f
356		!! Torsion array provided by C, dimensioned by getNlocal
357	<	real( kind = dp ), dimension(3,getNlocal()) :: t
357	>	real( kind = dp ), dimension(3,nLocal) :: t
358	>
359		!! Stress Tensor
360		real( kind = dp), dimension(9) :: tau
361		real ( kind = dp ) :: pot
#	Line 217 \| Line 368 \| contains
368		integer :: ncol
369		integer :: nprocs
370		#endif
220	–	integer :: nlocal
371		integer :: natoms
372		logical :: update_nlist
373		integer :: i, j, jbeg, jend, jnab
#	Line 225 \| Line 375 \| contains
375		real( kind = DP ) :: rijsq
376		real(kind=dp),dimension(3) :: d
377		real(kind=dp) :: rfpot, mu_i, virial
378	<	integer :: me_i
378	>	integer :: me_i, me_j
379		logical :: is_dp_i
380		integer :: neighborListSize
381		integer :: listerror, error
382		integer :: localError
383	+	integer :: propPack_i, propPack_j
384
385	<	real(kind=dp) :: listSkin = 1.0
235	<
385	>	real(kind=dp) :: listSkin = 1.0
386
387		!! initialize local variables
388
389		#ifdef IS_MPI
390		pot_local = 0.0_dp
241	–	nlocal = getNlocal()
391		nrow = getNrow(plan_row)
392		ncol = getNcol(plan_col)
393		#else
245	–	nlocal = getNlocal()
394		natoms = nlocal
395		#endif
396
397	<	call check_initialization(localError)
397	>	call doReadyCheck(localError)
398		if ( localError .ne. 0 ) then
399	<	call handleError("do_force_loop","Not Initialized")
399	>	call handleError("do_force_loop", "Not Initialized")
400		error = -1
401		return
402		end if
#	Line 257 \| Line 405 \| contains
405		do_pot = do_pot_c
406		do_stress = do_stress_c
407
260	–
408		! Gather all information needed by all force loops:
409
410		#ifdef IS_MPI
#	Line 265 \| Line 412 \| contains
412		call gather(q,q_Row,plan_row3d)
413		call gather(q,q_Col,plan_col3d)
414
415	<	if (FF_UsesDirectionalAtoms() .and. SimUsesDirectionalAtoms()) then
415	>	if (FF_UsesDirectionalAtoms() .and. SIM_uses_directional_atoms) then
416		call gather(u_l,u_l_Row,plan_row3d)
417		call gather(u_l,u_l_Col,plan_col3d)
418
#	Line 281 \| Line 428 \| contains
428		nloops = nloops + 1
429		#endif
430
431	<	if (FF_RequiresPrepairCalc() .and. SimRequiresPrepairCalc()) then
431	>	if (FF_RequiresPrepairCalc() .and. SIM_requires_prepair_calc) then
432		!! See if we need to update neighbor lists
433		call checkNeighborList(nlocal, q, listSkin, update_nlist)
434		!! if_mpi_gather_stuff_for_prepair
#	Line 445 \| Line 592 \| contains
592		#ifdef IS_MPI
593
594		if (update_nlist) then
448	–
595		!! save current configuration, construct neighbor list,
596		!! and calculate forces
597		call saveNeighborList(nlocal, q)
#	Line 454 \| Line 600 \| contains
600		nlist = 0
601
602		do i = 1, nrow
603	+
604		point(i) = nlist + 1
605
606		inner: do j = 1, ncol
#	Line 511 \| Line 658 \| contains
658		#else
659
660		if (update_nlist) then
661	<
661	>
662		! save current configuration, contruct neighbor list,
663		! and calculate forces
664		call saveNeighborList(natoms, q)
#	Line 602 \| Line 749 \| contains
749		f(1:3,i) = f(1:3,i) + f_temp(1:3,i)
750		end do
751
752	<	if (FF_UsesDirectionalAtoms() .and. SimUsesDirectionalAtoms()) then
752	>	if (FF_UsesDirectionalAtoms() .and. SIM_uses_directional_atoms) then
753		t_temp = 0.0_dp
754		call scatter(t_Row,t_temp,plan_row3d)
755		do i = 1,nlocal
#	Line 636 \| Line 783 \| contains
783		endif
784		#endif
785
786	<	if (FF_RequiresPostpairCalc() .and. SimRequiresPostpairCalc()) then
786	>	if (FF_RequiresPostpairCalc() .and. SIM_requires_postpair_calc) then
787
788	<	if (FF_uses_RF .and. SimUsesRF()) then
788	>	if (FF_uses_RF .and. SIM_uses_RF) then
789
790		#ifdef IS_MPI
791		call scatter(rf_Row,rf,plan_row3d)
#	Line 648 \| Line 795 \| contains
795		end do
796		#endif
797
798	<	do i = 1, getNlocal()
798	>	do i = 1, nLocal
799
800		rfpot = 0.0_DP
801		#ifdef IS_MPI
#	Line 656 \| Line 803 \| contains
803		#else
804		me_i = atid(i)
805		#endif
806	<	call getElementProperty(atypes, me_i, "is_DP", is_DP_i)
807	<	if ( is_DP_i ) then
808	<	call getElementProperty(atypes, me_i, "dipole_moment", mu_i)
806	>
807	>	if (PropertyMap(me_i)%is_DP) then
808	>
809	>	mu_i = PropertyMap(me_i)%dipole_moment
810	>
811		!! The reaction field needs to include a self contribution
812		!! to the field:
813	<	call accumulate_self_rf(i, mu_i, u_l)
813	>	call accumulate_self_rf(i, mu_i, u_l)
814		!! Get the reaction field contribution to the
815		!! potential and torques:
816		call reaction_field_final(i, mu_i, u_l, rfpot, t, do_pot)
#	Line 698 \| Line 847 \| contains
847		tau = tau_Temp
848		virial = virial_Temp
849		endif
850	<
850	>
851		#endif
852	<
853	<	#ifdef PROFILE
854	<	if (do_pot) then
706	<
707	<	#ifdef IS_MPI
708	<
709	<
710	<	call printCommTime()
711	<
712	<	call mpi_allreduce(forceTime,globalForceTime,1,MPI_DOUBLE_PRECISION, &
713	<	mpi_sum,mpi_comm_world,mpi_err)
714	<
715	<	call mpi_allreduce(forceTime,maxForceTime,1,MPI_DOUBLE_PRECISION, &
716	<	MPI_MAX,mpi_comm_world,mpi_err)
717	<
718	<	call mpi_comm_size( MPI_COMM_WORLD, nprocs,mpi_err)
719	<
720	<	if (getMyNode() == 0) then
721	<	write(,) "Total processor time spent in force calculations is: ", globalForceTime
722	<	write(,) "Total Time spent in force loop per processor is: ", globalforceTime/nprocs
723	<	write(,) "Maximum force time on any processor is: ", maxForceTime
724	<	end if
725	<	#else
726	<	write(,) "Time spent in force loop is: ", forceTime
727	<	#endif
728	<
729	<
730	<	endif
731	<
732	<	#endif
733	<
852	>
853	>
854	>
855		end subroutine do_force_loop
856
857		subroutine do_pair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
858
859		real( kind = dp ) :: pot
860	<	real( kind = dp ), dimension(3,getNlocal()) :: u_l
861	<	real (kind=dp), dimension(9,getNlocal()) :: A
862	<	real (kind=dp), dimension(3,getNlocal()) :: f
863	<	real (kind=dp), dimension(3,getNlocal()) :: t
860	>	real( kind = dp ), dimension(3,nLocal) :: u_l
861	>	real (kind=dp), dimension(9,nLocal) :: A
862	>	real (kind=dp), dimension(3,nLocal) :: f
863	>	real (kind=dp), dimension(3,nLocal) :: t
864
865		logical, intent(inout) :: do_pot, do_stress
866		integer, intent(in) :: i, j
#	Line 752 \| Line 873 \| contains
873		logical :: is_EAM_i,is_EAM_j
874		logical :: is_Sticky_i, is_Sticky_j
875		integer :: me_i, me_j
876	<
876	>	integer :: propPack_i
877	>	integer :: propPack_j
878		r = sqrt(rijsq)
879
880		#ifdef IS_MPI
#	Line 769 \| Line 891 \| contains
891		me_j = atid(j)
892
893		#endif
894	+
895	+	if (FF_uses_LJ .and. SIM_uses_LJ) then
896
897	<	if (FF_uses_LJ .and. SimUsesLJ()) then
774	<	call getElementProperty(atypes, me_i, "is_LJ", is_LJ_i)
775	<	call getElementProperty(atypes, me_j, "is_LJ", is_LJ_j)
776	<
777	<	if ( is_LJ_i .and. is_LJ_j ) &
897	>	if ( PropertyMap(me_i)%is_LJ .and. PropertyMap(me_j)%is_LJ ) &
898		call do_lj_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
899	+
900		endif
901
902	<	if (FF_uses_dipoles .and. SimUsesDipoles()) then
782	<	call getElementProperty(atypes, me_i, "is_DP", is_DP_i)
783	<	call getElementProperty(atypes, me_j, "is_DP", is_DP_j)
902	>	if (FF_uses_dipoles .and. SIM_uses_dipoles) then
903
904	<	if ( is_DP_i .and. is_DP_j ) then
904	>	if ( PropertyMap(me_i)%is_DP .and. PropertyMap(me_j)%is_DP) then
905		call do_dipole_pair(i, j, d, r, rijsq, pot, u_l, f, t, &
906		do_pot, do_stress)
907	<	if (FF_uses_RF .and. SimUsesRF()) then
907	>	if (FF_uses_RF .and. SIM_uses_RF) then
908		call accumulate_rf(i, j, r, u_l)
909		call rf_correct_forces(i, j, d, r, u_l, f, do_stress)
910		endif
#	Line 793 \| Line 912 \| contains
912		endif
913		endif
914
915	<	if (FF_uses_Sticky .and. SimUsesSticky()) then
915	>	if (FF_uses_Sticky .and. SIM_uses_sticky) then
916
917	<	call getElementProperty(atypes, me_i, "is_Sticky", is_Sticky_i)
799	<	call getElementProperty(atypes, me_j, "is_Sticky", is_Sticky_j)
800	<
801	<	if ( is_Sticky_i .and. is_Sticky_j ) then
917	>	if ( PropertyMap(me_i)%is_Sticky .and. PropertyMap(me_j)%is_Sticky) then
918		call do_sticky_pair(i, j, d, r, rijsq, A, pot, f, t, &
919		do_pot, do_stress)
920		endif
921		endif
922
923
924	<	if (FF_uses_GB .and. SimUsesGB()) then
809	<
810	<
811	<	call getElementProperty(atypes, me_i, "is_GB", is_GB_i)
812	<	call getElementProperty(atypes, me_j, "is_GB", is_GB_j)
924	>	if (FF_uses_GB .and. SIM_uses_GB) then
925
926	<	if ( is_GB_i .and. is_GB_j ) then
926	>	if ( PropertyMap(me_i)%is_GB .and. PropertyMap(me_j)%is_GB) then
927		call do_gb_pair(i, j, d, r, rijsq, u_l, pot, f, t, &
928		do_pot, do_stress)
929		endif
930	+
931		endif
932
933
934
935	<	if (FF_uses_EAM .and. SimUsesEAM()) then
823	<	call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
824	<	call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
935	>	if (FF_uses_EAM .and. SIM_uses_EAM) then
936
937	<	if ( is_EAM_i .and. is_EAM_j ) &
938	<	call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
937	>	if ( PropertyMap(me_i)%is_EAM .and. PropertyMap(me_j)%is_EAM) then
938	>	call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
939	>	endif
940	>
941		endif
829	–
942
831	–
832	–
943		end subroutine do_pair
944
945
946
947		subroutine do_prepair(i, j, rijsq, d, do_pot, do_stress, u_l, A, f, t, pot)
948		real( kind = dp ) :: pot
949	<	real( kind = dp ), dimension(3,getNlocal()) :: u_l
950	<	real (kind=dp), dimension(9,getNlocal()) :: A
951	<	real (kind=dp), dimension(3,getNlocal()) :: f
952	<	real (kind=dp), dimension(3,getNlocal()) :: t
949	>	real( kind = dp ), dimension(3,nLocal) :: u_l
950	>	real (kind=dp), dimension(9,nLocal) :: A
951	>	real (kind=dp), dimension(3,nLocal) :: f
952	>	real (kind=dp), dimension(3,nLocal) :: t
953
954		logical, intent(inout) :: do_pot, do_stress
955		integer, intent(in) :: i, j
#	Line 869 \| Line 979 \| contains
979
980		#endif
981
982	<	if (FF_uses_EAM .and. SimUsesEAM()) then
983	<	call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
984	<	call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
875	<
876	<	if ( is_EAM_i .and. is_EAM_j ) &
982	>	if (FF_uses_EAM .and. SIM_uses_EAM) then
983	>
984	>	if (PropertyMap(me_i)%is_EAM .and. PropertyMap(me_j)%is_EAM) &
985		call calc_EAM_prepair_rho(i, j, d, r, rijsq )
878	–	endif
986
987	+	endif
988	+
989		end subroutine do_prepair
990
991
#	Line 886 \| Line 995 \| contains
995		integer :: nlocal
996		real( kind = dp ) :: pot
997
998	<	if (FF_uses_EAM .and. SimUsesEAM()) then
998	>	if (FF_uses_EAM .and. SIM_uses_EAM) then
999		call calc_EAM_preforce_Frho(nlocal,pot)
1000		endif
1001
#	Line 905 \| Line 1014 \| contains
1014		d(1:3) = q_j(1:3) - q_i(1:3)
1015
1016		! Wrap back into periodic box if necessary
1017	<	if ( SimUsesPBC() ) then
1017	>	if ( SIM_uses_PBC ) then
1018
1019		if( .not.boxIsOrthorhombic ) then
1020		! calc the scaled coordinates.
#	Line 945 \| Line 1054 \| contains
1054
1055		end subroutine get_interatomic_vector
1056
948	–	subroutine check_initialization(error)
949	–	integer, intent(out) :: error
950	–
951	–	error = 0
952	–	! Make sure we are properly initialized.
953	–	if (.not. do_forces_initialized) then
954	–	write(,) "Forces not initialized"
955	–	error = -1
956	–	return
957	–	endif
958	–
959	–	#ifdef IS_MPI
960	–	if (.not. isMPISimSet()) then
961	–	write(default_error,*) "ERROR: mpiSimulation has not been initialized!"
962	–	error = -1
963	–	return
964	–	endif
965	–	#endif
966	–
967	–	return
968	–	end subroutine check_initialization
969	–
970	–
1057		subroutine zero_work_arrays()
1058
1059		#ifdef IS_MPI
#	Line 1000 \| Line 1086 \| contains
1086		#endif
1087
1088
1089	<	if (FF_uses_EAM .and. SimUsesEAM()) then
1089	>	if (FF_uses_EAM .and. SIM_uses_EAM) then
1090		call clean_EAM()
1091		endif
1092
#	Line 1120 \| Line 1206 \| contains
1206		doesit = FF_uses_RF
1207		end function FF_RequiresPostpairCalc
1208
1209	+	#ifdef PROFILE
1210	+	function getforcetime() result(totalforcetime)
1211	+	real(kind=dp) :: totalforcetime
1212	+	totalforcetime = forcetime
1213	+	end function getforcetime
1214	+	#endif
1215	+
1216		!! This cleans componets of force arrays belonging only to fortran
1217
1218		end module do_Forces

Diff Legend

-–
+Removed lines
-+
+Added lines
-<
+Changed lines
->
+Changed lines

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents): Revision 843 by mmeineke, Wed Oct 29 20:41:39 2003 UTC vs. Revision 900 by chuckv, Tue Jan 6 18:54:57 2004 UTC

Diff Legend

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents):
Revision 843 by mmeineke, Wed Oct 29 20:41:39 2003 UTC vs.
Revision 900 by chuckv, Tue Jan 6 18:54:57 2004 UTC