ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/do_Forces.F90
(Generate patch)

Comparing trunk/OOPSE/libmdtools/do_Forces.F90 (file contents):
Revision 898 by chuckv, Mon Jan 5 22:49:14 2004 UTC vs.
Revision 946 by gezelter, Thu Jan 15 04:33:24 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.44 2004-01-05 22:49:14 chuckv Exp $, $Date: 2004-01-05 22:49:14 $, $Name: not supported by cvs2svn $, $Revision: 1.44 $
7 > !! @version $Id: do_Forces.F90,v 1.47 2004-01-15 04:33:24 gezelter Exp $, $Date: 2004-01-15 04:33:24 $, $Name: not supported by cvs2svn $, $Revision: 1.47 $
8  
9   module do_Forces
10    use force_globals
# Line 15 | Line 15 | module do_Forces
15    use lj
16    use sticky_pair
17    use dipole_dipole
18 +  use charge_charge
19    use reaction_field
20    use gb_pair
21    use vector_class
# Line 30 | Line 31 | module do_Forces
31   #define __FORTRAN90
32   #include "fForceField.h"
33  
34 <  logical, save :: do_forces_initialized = .false., haveRlist = .false.
34 >  logical, save :: haveRlist = .false.
35 >  logical, save :: haveNeighborList = .false.
36    logical, save :: havePolicies = .false.
37 +  logical, save :: haveSIMvariables = .false.
38 +  logical, save :: havePropertyMap = .false.
39 +  logical, save :: haveSaneForceField = .false.
40    logical, save :: FF_uses_LJ
41    logical, save :: FF_uses_sticky
42 +  logical, save :: FF_uses_charges
43    logical, save :: FF_uses_dipoles
44    logical, save :: FF_uses_RF
45    logical, save :: FF_uses_GB
46    logical, save :: FF_uses_EAM
47 +  logical, save :: SIM_uses_LJ
48 +  logical, save :: SIM_uses_sticky
49 +  logical, save :: SIM_uses_charges
50 +  logical, save :: SIM_uses_dipoles
51 +  logical, save :: SIM_uses_RF
52 +  logical, save :: SIM_uses_GB
53 +  logical, save :: SIM_uses_EAM
54 +  logical, save :: SIM_requires_postpair_calc
55 +  logical, save :: SIM_requires_prepair_calc
56 +  logical, save :: SIM_uses_directional_atoms
57 +  logical, save :: SIM_uses_PBC
58  
59    real(kind=dp), save :: rlist, rlistsq
60  
# Line 45 | Line 62 | module do_Forces
62    public :: do_force_loop
63    public :: setRlistDF
64  
65 +
66   #ifdef PROFILE
67    public :: getforcetime
68    real, save :: forceTime = 0
# Line 52 | Line 70 | module do_Forces
70    integer :: nLoops
71   #endif
72  
73 <  logical, allocatable :: propertyMapI(:,:)
74 <  logical, allocatable :: propertyMapJ(:,:)
73 >  type :: Properties
74 >     logical :: is_lj     = .false.
75 >     logical :: is_sticky = .false.
76 >     logical :: is_dp     = .false.
77 >     logical :: is_gb     = .false.
78 >     logical :: is_eam    = .false.
79 >     logical :: is_charge = .false.
80 >     real(kind=DP) :: charge = 0.0_DP
81 >     real(kind=DP) :: dipole_moment = 0.0_DP
82 >  end type Properties
83  
84 +  type(Properties), dimension(:),allocatable :: PropertyMap
85 +
86   contains
87  
88    subroutine setRlistDF( this_rlist )
# Line 65 | Line 93 | contains
93      rlistsq = rlist * rlist
94      
95      haveRlist = .true.
68    if( havePolicies ) do_forces_initialized = .true.
96  
97    end subroutine setRlistDF    
98 +
99 +  subroutine createPropertyMap(status)
100 +    integer :: nAtypes
101 +    integer :: status
102 +    integer :: i
103 +    logical :: thisProperty
104 +    real (kind=DP) :: thisDPproperty
105 +
106 +    status = 0
107 +
108 +    nAtypes = getSize(atypes)
109 +
110 +    if (nAtypes == 0) then
111 +       status = -1
112 +       return
113 +    end if
114 +        
115 +    if (.not. allocated(PropertyMap)) then
116 +       allocate(PropertyMap(nAtypes))
117 +    endif
118 +
119 +    do i = 1, nAtypes
120 +       call getElementProperty(atypes, i, "is_LJ", thisProperty)
121 +       PropertyMap(i)%is_LJ = thisProperty
122 +
123 +       call getElementProperty(atypes, i, "is_Charge", thisProperty)
124 +       PropertyMap(i)%is_Charge = thisProperty
125 +      
126 +       if (thisProperty) then
127 +          call getElementProperty(atypes, i, "charge", thisDPproperty)
128 +          PropertyMap(i)%charge = thisDPproperty
129 +       endif
130 +
131 +       call getElementProperty(atypes, i, "is_DP", thisProperty)
132 +       PropertyMap(i)%is_DP = thisProperty
133 +
134 +       if (thisProperty) then
135 +          call getElementProperty(atypes, i, "dipole_moment", thisDPproperty)
136 +          PropertyMap(i)%dipole_moment = thisDPproperty
137 +       endif
138 +
139 +       call getElementProperty(atypes, i, "is_Sticky", thisProperty)
140 +       PropertyMap(i)%is_Sticky = thisProperty
141 +       call getElementProperty(atypes, i, "is_GB", thisProperty)
142 +       PropertyMap(i)%is_GB = thisProperty
143 +       call getElementProperty(atypes, i, "is_EAM", thisProperty)
144 +       PropertyMap(i)%is_EAM = thisProperty
145 +    end do
146 +
147 +    havePropertyMap = .true.
148 +
149 +  end subroutine createPropertyMap
150 +
151 +  subroutine setSimVariables()
152 +    SIM_uses_LJ = SimUsesLJ()
153 +    SIM_uses_sticky = SimUsesSticky()
154 +    SIM_uses_charges = SimUsesCharges()
155 +    SIM_uses_dipoles = SimUsesDipoles()
156 +    SIM_uses_RF = SimUsesRF()
157 +    SIM_uses_GB = SimUsesGB()
158 +    SIM_uses_EAM = SimUsesEAM()
159 +    SIM_requires_postpair_calc = SimRequiresPostpairCalc()
160 +    SIM_requires_prepair_calc = SimRequiresPrepairCalc()
161 +    SIM_uses_directional_atoms = SimUsesDirectionalAtoms()
162 +    SIM_uses_PBC = SimUsesPBC()
163 +
164 +    haveSIMvariables = .true.
165 +
166 +    return
167 +  end subroutine setSimVariables
168  
169 +  subroutine doReadyCheck(error)
170 +    integer, intent(out) :: error
171 +
172 +    integer :: myStatus
173 +
174 +    error = 0
175 +    
176 +    if (.not. havePropertyMap) then
177 +
178 +       myStatus = 0
179 +
180 +       call createPropertyMap(myStatus)
181 +
182 +       if (myStatus .ne. 0) then
183 +          write(default_error, *) 'createPropertyMap failed in do_Forces!'
184 +          error = -1
185 +          return
186 +       endif
187 +    endif
188 +
189 +    if (.not. haveSIMvariables) then
190 +       call setSimVariables()
191 +    endif
192 +
193 +    if (.not. haveRlist) then
194 +       write(default_error, *) 'rList has not been set in do_Forces!'
195 +       error = -1
196 +       return
197 +    endif
198 +
199 +    if (SIM_uses_LJ .and. FF_uses_LJ) then
200 +       if (.not. havePolicies) then
201 +          write(default_error, *) 'LJ mixing Policies have not been set in do_Forces!'
202 +          error = -1
203 +          return
204 +       endif
205 +    endif
206 +
207 +    if (.not. haveNeighborList) then
208 +       write(default_error, *) 'neighbor list has not been initialized in do_Forces!'
209 +       error = -1
210 +       return
211 +    end if
212 +
213 +    if (.not. haveSaneForceField) then
214 +       write(default_error, *) 'Force Field is not sane in do_Forces!'
215 +       error = -1
216 +       return
217 +    end if
218 +
219 + #ifdef IS_MPI
220 +    if (.not. isMPISimSet()) then
221 +       write(default_error,*) "ERROR: mpiSimulation has not been initialized!"
222 +       error = -1
223 +       return
224 +    endif
225 + #endif
226 +    return
227 +  end subroutine doReadyCheck
228 +    
229 +
230    subroutine init_FF(LJMIXPOLICY, use_RF_c, thisStat)
231  
232      integer, intent(in) :: LJMIXPOLICY
# Line 93 | Line 251 | contains
251    
252      FF_uses_LJ = .false.
253      FF_uses_sticky = .false.
254 +    FF_uses_charges = .false.
255      FF_uses_dipoles = .false.
256      FF_uses_GB = .false.
257      FF_uses_EAM = .false.
258      
259      call getMatchingElementList(atypes, "is_LJ", .true., nMatches, MatchList)
260      if (nMatches .gt. 0) FF_uses_LJ = .true.
261 +
262 +    call getMatchingElementList(atypes, "is_Charge", .true., nMatches, MatchList)
263 +    if (nMatches .gt. 0) FF_uses_charges = .true.
264      
265      call getMatchingElementList(atypes, "is_DP", .true., nMatches, MatchList)
266      if (nMatches .gt. 0) FF_uses_dipoles = .true.
# Line 113 | Line 275 | contains
275      call getMatchingElementList(atypes, "is_EAM", .true., nMatches, MatchList)
276      if (nMatches .gt. 0) FF_uses_EAM = .true.
277      
278 +    !! Assume sanity (for the sake of argument)
279 +    haveSaneForceField = .true.
280 +
281      !! check to make sure the FF_uses_RF setting makes sense
282      
283      if (FF_uses_dipoles) then
# Line 124 | Line 289 | contains
289         if (FF_uses_RF) then          
290            write(default_error,*) 'Using Reaction Field with no dipoles?  Huh?'
291            thisStat = -1
292 +          haveSaneForceField = .false.
293            return
294         endif
295      endif
# Line 138 | Line 304 | contains
304         case default
305            write(default_error,*) 'unknown LJ Mixing Policy!'
306            thisStat = -1
307 +          haveSaneForceField = .false.
308            return            
309         end select
310         if (my_status /= 0) then
311            thisStat = -1
312 +          haveSaneForceField = .false.
313            return
314         end if
315 +       havePolicies = .true.
316      endif
317  
318      if (FF_uses_sticky) then
319         call check_sticky_FF(my_status)
320         if (my_status /= 0) then
321            thisStat = -1
322 +          haveSaneForceField = .false.
323            return
324         end if
325      endif
# Line 158 | Line 328 | contains
328      if (FF_uses_EAM) then
329           call init_EAM_FF(my_status)
330         if (my_status /= 0) then
331 <          write(*,*) "init_EAM_FF returned a bad status"
331 >          write(default_error, *) "init_EAM_FF returned a bad status"
332            thisStat = -1
333 +          haveSaneForceField = .false.
334            return
335         end if
336      endif
337  
167
168    
338      if (FF_uses_GB) then
339         call check_gb_pair_FF(my_status)
340         if (my_status .ne. 0) then
341            thisStat = -1
342 +          haveSaneForceField = .false.
343            return
344         endif
345      endif
346  
347      if (FF_uses_GB .and. FF_uses_LJ) then
348      endif
349 <    if (.not. do_forces_initialized) then
349 >    if (.not. haveNeighborList) then
350         !! Create neighbor lists
351         call expandNeighborList(nLocal, my_status)
352         if (my_Status /= 0) then
# Line 184 | Line 354 | contains
354            thisStat = -1
355            return
356         endif
357 +       haveNeighborList = .true.
358      endif
359      
189
190    havePolicies = .true.
191    if( haveRlist ) do_forces_initialized = .true.
192
360    end subroutine init_FF
361    
362  
# Line 246 | Line 413 | contains
413      natoms = nlocal
414   #endif
415  
416 <    call check_initialization(localError)
416 >    call doReadyCheck(localError)
417      if ( localError .ne. 0 ) then
418 <       call handleError("do_force_loop","Not Initialized")
418 >       call handleError("do_force_loop", "Not Initialized")
419         error = -1
420         return
421      end if
# Line 256 | Line 423 | contains
423  
424      do_pot = do_pot_c
425      do_stress = do_stress_c
259
260
261 #ifdef IS_MPI
262    if (.not.allocated(propertyMapI)) then
263       allocate(propertyMapI(5,nrow))
264    endif
265
266    do i = 1, nrow
267       me_i = atid_row(i)
268 #else
269    if (.not.allocated(propertyMapI)) then
270       allocate(propertyMapI(5,nlocal))
271    endif
426  
273    do i = 1, natoms
274       me_i = atid(i)
275 #endif
276      
277       propertyMapI(1:5,i) = .false.
278
279       call getElementProperty(atypes, me_i, "propertyPack", propPack_i)
280    
281       ! unpack the properties
282      
283       if (iand(propPack_i, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
284            propertyMapI(1, i) = .true.
285       if (iand(propPack_i, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
286            propertyMapI(2, i) = .true.
287       if (iand(propPack_i, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
288            propertyMapI(3, i) = .true.
289       if (iand(propPack_i, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
290            propertyMapI(4, i) = .true.
291       if (iand(propPack_i, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
292            propertyMapI(5, i) = .true.
293
294    end do
295
296 #ifdef IS_MPI
297    if (.not.allocated(propertyMapJ)) then
298       allocate(propertyMapJ(5,ncol))
299    endif
300
301    do j = 1, ncol
302       me_j = atid_col(j)
303 #else
304    if (.not.allocated(propertyMapJ)) then
305       allocate(propertyMapJ(5,nlocal))
306    endif
307
308    do j = 1, natoms
309       me_j = atid(j)
310 #endif
311      
312       propertyMapJ(1:5,j) = .false.
313
314       call getElementProperty(atypes, me_j, "propertyPack", propPack_j)
315    
316       ! unpack the properties
317      
318       if (iand(propPack_j, LJ_PROPERTY_MASK) .eq. LJ_PROPERTY_MASK) &
319            propertyMapJ(1, j) = .true.
320       if (iand(propPack_j, DP_PROPERTY_MASK) .eq. DP_PROPERTY_MASK) &
321            propertyMapJ(2, j) = .true.
322       if (iand(propPack_j, STICKY_PROPERTY_MASK) .eq. STICKY_PROPERTY_MASK) &
323            propertyMapJ(3, j) = .true.
324       if (iand(propPack_j, GB_PROPERTY_MASK) .eq. GB_PROPERTY_MASK) &
325            propertyMapJ(4, j) = .true.
326       if (iand(propPack_j, EAM_PROPERTY_MASK) .eq. EAM_PROPERTY_MASK) &
327            propertyMapJ(5, j) = .true.
328
329    end do
330
427      ! Gather all information needed by all force loops:
428      
429   #ifdef IS_MPI    
# Line 335 | Line 431 | contains
431      call gather(q,q_Row,plan_row3d)
432      call gather(q,q_Col,plan_col3d)
433          
434 <    if (FF_UsesDirectionalAtoms() .and. SimUsesDirectionalAtoms()) then
434 >    if (FF_UsesDirectionalAtoms() .and. SIM_uses_directional_atoms) then
435         call gather(u_l,u_l_Row,plan_row3d)
436         call gather(u_l,u_l_Col,plan_col3d)
437        
# Line 351 | Line 447 | contains
447      nloops = nloops + 1
448   #endif
449    
450 <    if (FF_RequiresPrepairCalc() .and. SimRequiresPrepairCalc()) then
450 >    if (FF_RequiresPrepairCalc() .and. SIM_requires_prepair_calc) then
451         !! See if we need to update neighbor lists
452         call checkNeighborList(nlocal, q, listSkin, update_nlist)  
453         !! if_mpi_gather_stuff_for_prepair
# Line 672 | Line 768 | contains
768         f(1:3,i) = f(1:3,i) + f_temp(1:3,i)
769      end do
770      
771 <    if (FF_UsesDirectionalAtoms() .and. SimUsesDirectionalAtoms()) then
771 >    if (FF_UsesDirectionalAtoms() .and. SIM_uses_directional_atoms) then
772         t_temp = 0.0_dp
773         call scatter(t_Row,t_temp,plan_row3d)
774         do i = 1,nlocal
# Line 706 | Line 802 | contains
802      endif    
803   #endif
804  
805 <    if (FF_RequiresPostpairCalc() .and. SimRequiresPostpairCalc()) then
805 >    if (FF_RequiresPostpairCalc() .and. SIM_requires_postpair_calc) then
806        
807 <       if (FF_uses_RF .and. SimUsesRF()) then
807 >       if (FF_uses_RF .and. SIM_uses_RF) then
808            
809   #ifdef IS_MPI
810            call scatter(rf_Row,rf,plan_row3d)
# Line 726 | Line 822 | contains
822   #else
823               me_i = atid(i)
824   #endif
825 <             call getElementProperty(atypes, me_i, "is_DP", is_DP_i)      
826 <             if ( is_DP_i ) then
827 <                call getElementProperty(atypes, me_i, "dipole_moment", mu_i)
825 >
826 >             if (PropertyMap(me_i)%is_DP) then
827 >
828 >                mu_i = PropertyMap(me_i)%dipole_moment
829 >
830                  !! The reaction field needs to include a self contribution
831                  !! to the field:
832 <                call accumulate_self_rf(i, mu_i, u_l)            
832 >                call accumulate_self_rf(i, mu_i, u_l)
833                  !! Get the reaction field contribution to the
834                  !! potential and torques:
835                  call reaction_field_final(i, mu_i, u_l, rfpot, t, do_pot)
# Line 788 | Line 886 | contains
886      real ( kind = dp ), intent(inout)    :: rijsq
887      real ( kind = dp )                :: r
888      real ( kind = dp ), intent(inout) :: d(3)
791    logical :: is_LJ_i, is_LJ_j
792    logical :: is_DP_i, is_DP_j
793    logical :: is_GB_i, is_GB_j
794    logical :: is_EAM_i,is_EAM_j
795    logical :: is_Sticky_i, is_Sticky_j
889      integer :: me_i, me_j
890 <    integer :: propPack_i
798 <    integer :: propPack_j
890 >
891      r = sqrt(rijsq)
892  
893   #ifdef IS_MPI
894      if (tagRow(i) .eq. tagColumn(j)) then
895         write(0,*) 'do_pair is doing', i , j, tagRow(i), tagColumn(j)
896      endif
805
897      me_i = atid_row(i)
898      me_j = atid_col(j)
808
899   #else
810
900      me_i = atid(i)
901      me_j = atid(j)
813
902   #endif
903      
904 <    if (FF_uses_LJ .and. SimUsesLJ()) then
905 <
906 <       if ( propertyMapI(1, me_i) .and. propertyMapJ(1, me_j) ) &
907 <            call do_lj_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
908 <
904 >    if (FF_uses_LJ .and. SIM_uses_LJ) then
905 >      
906 >       if ( PropertyMap(me_i)%is_LJ .and. PropertyMap(me_j)%is_LJ ) then
907 >          call do_lj_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
908 >       endif
909 >      
910      endif
911 <
912 <    if (FF_uses_dipoles .and. SimUsesDipoles()) then
911 >    
912 >    if (FF_uses_charges .and. SIM_uses_charges) then
913        
914 <       if ( propertyMapI(2, me_i) .and. propertyMapJ(2, me_j)) then
914 >       if (PropertyMap(me_i)%is_Charge .and. PropertyMap(me_j)%is_Charge) then
915 >          call do_charge_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
916 >       endif
917 >      
918 >    endif
919 >    
920 >    if (FF_uses_dipoles .and. SIM_uses_dipoles) then
921 >      
922 >       if ( PropertyMap(me_i)%is_DP .and. PropertyMap(me_j)%is_DP) then
923            call do_dipole_pair(i, j, d, r, rijsq, pot, u_l, f, t, &
924                 do_pot, do_stress)
925 <          if (FF_uses_RF .and. SimUsesRF()) then
925 >          if (FF_uses_RF .and. SIM_uses_RF) then
926               call accumulate_rf(i, j, r, u_l)
927               call rf_correct_forces(i, j, d, r, u_l, f, do_stress)
928 <          endif
832 <          
928 >          endif          
929         endif
930 +
931      endif
932  
933 <    if (FF_uses_Sticky .and. SimUsesSticky()) then
933 >    if (FF_uses_Sticky .and. SIM_uses_sticky) then
934  
935 <       if ( propertyMapI(3, me_i) .and. propertyMapJ(3, me_j)) then
935 >       if ( PropertyMap(me_i)%is_Sticky .and. PropertyMap(me_j)%is_Sticky) then
936            call do_sticky_pair(i, j, d, r, rijsq, A, pot, f, t, &
937                 do_pot, do_stress)
938         endif
939 +
940      endif
941  
942  
943 <    if (FF_uses_GB .and. SimUsesGB()) then
943 >    if (FF_uses_GB .and. SIM_uses_GB) then
944        
945 <       if ( propertyMapI(4, me_i) .and. propertyMapJ(4, me_j)) then
945 >       if ( PropertyMap(me_i)%is_GB .and. PropertyMap(me_j)%is_GB) then
946            call do_gb_pair(i, j, d, r, rijsq, u_l, pot, f, t, &
947                 do_pot, do_stress)          
948         endif
949  
950      endif
951 <    
952 <
953 <  
954 <   if (FF_uses_EAM .and. SimUsesEAM()) then
955 <      
956 <      if ( propertyMapI(5, me_i) .and. propertyMapJ(5, me_j)) then
957 <         call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
958 <      endif
861 <
862 <   endif
951 >      
952 >    if (FF_uses_EAM .and. SIM_uses_EAM) then
953 >      
954 >       if ( PropertyMap(me_i)%is_EAM .and. PropertyMap(me_j)%is_EAM) then
955 >          call do_eam_pair(i, j, d, r, rijsq, pot, f, do_pot, do_stress)
956 >       endif
957 >      
958 >    endif
959  
960    end subroutine do_pair
961  
# Line 900 | Line 996 | contains
996    
997   #endif
998      
999 <   if (FF_uses_EAM .and. SimUsesEAM()) then
1000 <      call getElementProperty(atypes, me_i, "is_EAM", is_EAM_i)
1001 <      call getElementProperty(atypes, me_j, "is_EAM", is_EAM_j)
906 <      
907 <      if ( is_EAM_i .and. is_EAM_j ) &
999 >   if (FF_uses_EAM .and. SIM_uses_EAM) then
1000 >
1001 >      if (PropertyMap(me_i)%is_EAM .and. PropertyMap(me_j)%is_EAM) &
1002             call calc_EAM_prepair_rho(i, j, d, r, rijsq )
909   endif
1003  
1004 +   endif
1005 +  
1006   end subroutine do_prepair
1007  
1008  
# Line 917 | Line 1012 | contains
1012      integer :: nlocal
1013      real( kind = dp ) :: pot
1014  
1015 <    if (FF_uses_EAM .and. SimUsesEAM()) then
1015 >    if (FF_uses_EAM .and. SIM_uses_EAM) then
1016         call calc_EAM_preforce_Frho(nlocal,pot)
1017      endif
1018  
# Line 936 | Line 1031 | contains
1031      d(1:3) = q_j(1:3) - q_i(1:3)
1032  
1033      ! Wrap back into periodic box if necessary
1034 <    if ( SimUsesPBC() ) then
1034 >    if ( SIM_uses_PBC ) then
1035        
1036         if( .not.boxIsOrthorhombic ) then
1037            ! calc the scaled coordinates.
# Line 976 | Line 1071 | contains
1071      
1072    end subroutine get_interatomic_vector
1073    
979  subroutine check_initialization(error)
980    integer, intent(out) :: error
981    
982    error = 0
983    ! Make sure we are properly initialized.
984    if (.not. do_forces_initialized) then
985       write(*,*) "Forces not initialized"
986       error = -1
987       return
988    endif
989
990 #ifdef IS_MPI
991    if (.not. isMPISimSet()) then
992       write(default_error,*) "ERROR: mpiSimulation has not been initialized!"
993       error = -1
994       return
995    endif
996 #endif
997    
998    return
999  end subroutine check_initialization
1000
1001  
1074    subroutine zero_work_arrays()
1075      
1076   #ifdef IS_MPI
# Line 1031 | Line 1103 | contains
1103   #endif
1104  
1105  
1106 <    if (FF_uses_EAM .and. SimUsesEAM()) then
1106 >    if (FF_uses_EAM .and. SIM_uses_EAM) then
1107         call clean_EAM()
1108      endif
1109  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines