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

Comparing trunk/OOPSE/libmdtools/f_verlet_constrained.F90 (file contents):
Revision 389 by mmeineke, Mon Mar 24 15:26:05 2003 UTC vs.
Revision 542 by mmeineke, Fri May 30 21:31:48 2003 UTC

# Line 1 | Line 1
1   subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
2       vx,vy,vz,fx,fy,fz, n_constrained, constraints_sqr, c_i, c_j, &
3 <     box_x, box_y, box_z)
3 >     box_x, box_y, box_z, isError)
4    implicit none
5    
6    ! *******************************************************************
# Line 16 | Line 16 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
16        
17    ! move part a calculate velocities
18    
19 <  INTEGER     I
19 >  INTEGER     I, isError
20    double precision DT2, DTSQ2
21    
22    double precision box_x, box_y, box_z
# Line 49 | Line 49 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
49    double precision pxab, pyab, pzab, pabsq, rabsq, diffsq
50    double precision rxab, ryab, rzab, rpab, gab
51    double precision dx, dy, dz, rma, rmb
52 +  double precision prSqr, rpabSqr;
53    integer a, b, it, maxit
54    parameter (rptol = 1.0d-6, tol = 1.0d-6 )
55    parameter ( maxit = 100 )
# Line 58 | Line 59 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
59    double precision e_convert
60    parameter ( e_convert = 4.184d-4 )
61    
62 +  isError = 0
63  
64    DT2   = DT / 2.0d0
65    DTSQ2 = DT * DT2
# Line 127 | Line 129 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
129  
130                rpab = rxab * pxab + ryab * pyab + rzab * pzab
131                
132 <              if( dabs(rpab) .lt. ( rabsq * rptol ) ) then
132 >              rpabSqr = rpab * rpab
133 >              prSqr = - diffsq
134 >
135 >              if( rpabSqr .lt. ( rabsq * prSqr )) then
136                   write (*, '('' Constraint Failure '')' )
137 <                 write (*,*) a-1, b-1,rpab, rabsq * rptol
138 <                 stop
137 >                 write (*,*) a-1, b-1,rpabSqr, rabsq * prSqr
138 >                 isError = 1
139 >                 return
140                end if
141                
142                rma = 1.0d0 / mass(a)
# Line 182 | Line 188 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
188  
189    if( .not. done) then
190       write (*, '('' too many constraint iterations in move_a '')' )
191 <     stop
191 >     isError = 1
192 >     return
193    endif
194  
195    ! store new values
# Line 194 | Line 201 | subroutine v_constrain_a(dt,natoms,mass,rx,ry,rz, &
201    enddo
202    
203    
197 –  
204    RETURN
205   end subroutine v_constrain_a
206  
207   Subroutine v_constrain_b(dt,natoms,mass,rx,ry,rz, &
208       vx,vy,vz,fx,fy,fz,k, n_constrained, constraints_sqr, &
209 <     c_i, c_j, box_x, box_y, box_z)
209 >     c_i, c_j, box_x, box_y, box_z, isError)
210    implicit none
211    
212    ! *******************************************************************
# Line 215 | Line 221 | Subroutine v_constrain_b(dt,natoms,mass,rx,ry,rz, &
221    
222    ! declarations
223    
224 <  integer  i
224 >  integer  i, isError
225    double precision accvel2, dt2
226    double precision box_x, box_y, box_z
227    
# Line 258 | Line 264 | Subroutine v_constrain_b(dt,natoms,mass,rx,ry,rz, &
264  
265    ! *******************************************************************
266    
267 +  isError = 0
268 +
269    tol = 1.0d-6 / dt
270    dt2 = dt / 2.0d0
271    accvel2 = 0.0d0
# Line 349 | Line 357 | Subroutine v_constrain_b(dt,natoms,mass,rx,ry,rz, &
357    if (.not. done) then
358      
359       write(*, '('' Too many constraint iterations in moveb '')')
360 <     stop
360 >     isError = 1
361 >     return
362      
363    endif
364    

Diff Legend

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