hdiff output

r31389/ljadd.f 2017-01-21 10:35:33.350250000 +0000 r31388/ljadd.f 2017-01-21 10:35:33.970250000 +0000
212: C212: C
213: C  Symmetrise Hessian213: C  Symmetrise Hessian
214: C214: C
215:       DO J1=1,3*N215:       DO J1=1,3*N
216:          DO J2=J1+1,3*N216:          DO J2=J1+1,3*N
217:             HESS(J1,J2)=HESS(J2,J1)217:             HESS(J1,J2)=HESS(J2,J1)
218:          ENDDO218:          ENDDO
219:       ENDDO219:       ENDDO
220:       RETURN220:       RETURN
221:       END221:       END
222: C 
223: C************************************************************************* 
224: C 
225: C  Subroutine LJADD2 calculates the cartesian gradient and second 
226: C  derivative matrix analytically for LJ with addressable epsilon values. Reduced units. 
227: C  This routine treats multiple copies of a target of cluster size NADDTARGET. 
228: C  The epsilon values are replicated via the MOD function. 
229: C 
230: C************************************************************************* 
231: C 
232:       SUBROUTINE LJADD2(N, X, V, ENERGY, GTEST, STEST) 
233:       USE COMMONS, ONLY : LJADDEPS, VT, NADDTARGET, MYUNIT 
234:       IMPLICIT NONE 
235:       INTEGER N, J1, J2, J3, J4, MJ1, MJ2 
236:       LOGICAL GTEST, STEST 
237:       DOUBLE PRECISION X(3*N), ENERGY, R6, 
238:      1                 V(3*N), R2(N,N), R2T, 
239:      2                 R8(N,N), G(N,N), XG(N,N), 
240:      3                 R14(N,N), F(N,N), DUMMY, DUMMYX, DUMMYY, DUMMYZ, DIST, XMUL2 
241: C  
242: C  Store distance matrices. 
243: C 
244: !     WRITE(MYUNIT,'(A)') 'coords in LJADD2:' 
245: !     WRITE(MYUNIT,'(3G20.10)') X(1:3*N) 
246:       ENERGY=0.0D0 
247:       DO J1=1,N 
248:          VT(J1)=0.0D0 
249:       ENDDO 
250:       IF (GTEST.AND.(.NOT.STEST)) THEN 
251:          DO J1=1,N 
252:             MJ1=MOD(J1,NADDTARGET)+1 
253:             J3=3*J1 
254:             XG(J1,J1)=0.0D0 
255:             DO J2=J1+1,N 
256:                MJ2=MOD(J2,NADDTARGET)+1 
257:                J4=3*J2 
258:                DIST=(X(J3-2)-X(J4-2))**2+(X(J3-1)-X(J4-1))**2+(X(J3)-X(J4))**2 
259:                DIST=1.0D0/DIST 
260:                R6=DIST**3 
261:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
262:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6 
263:                ELSE 
264:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0) 
265:                ENDIF 
266:                VT(J1)=VT(J1)+DUMMY 
267:                VT(J2)=VT(J2)+DUMMY 
268:                ENERGY=ENERGY+DUMMY 
269: !              WRITE(MYUNIT,'(A,2I8,2G20.10)') 'J1,J2,DUMMY,ENERGY=',J1,J2,DUMMY,ENERGY 
270:                DIST=DIST*R6 
271:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
272:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*DIST 
273:                ELSE 
274:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*DIST 
275:                ENDIF 
276:                XG(J1,J2)=XG(J2,J1) 
277:             ENDDO 
278:          ENDDO 
279:       ELSEIF (GTEST) THEN 
280:          DO J1=1,N 
281:             MJ1=MOD(J1,NADDTARGET)+1 
282:             XG(J1,J1)=0.0D0 
283:             R2(J1,J1)=0.0D0 
284:             R8(J1,J1)=0.0D0 
285:             R14(J1,J1)=0.0D0 
286:             DO J2=J1+1,N 
287:                MJ2=MOD(J2,NADDTARGET)+1 
288:                R2(J2,J1)=(X(3*(J1-1)+1)-X(3*(J2-1)+1))**2 
289:      1                  +(X(3*(J1-1)+2)-X(3*(J2-1)+2))**2 
290:      2                  +(X(3*(J1-1)+3)-X(3*(J2-1)+3))**2 
291:                R2(J2,J1)=1.0D0/R2(J2,J1) 
292:                R6=R2(J2,J1)**3 
293:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
294:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6 
295:                ELSE 
296:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0) 
297:                ENDIF 
298:                VT(J1)=VT(J1)+DUMMY 
299:                VT(J2)=VT(J2)+DUMMY 
300:                ENERGY=ENERGY+DUMMY 
301:                R8(J2,J1)=R2(J2,J1)**4 
302:                R14(J2,J1)=R8(J2,J1)*R8(J2,J1)/R2(J2,J1) 
303:                R2(J1,J2)=R2(J2,J1) 
304:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
305:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*R2(J1,J2)*R6 
306:                ELSE 
307:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*R2(J1,J2)*R6 
308:                ENDIF 
309:                XG(J1,J2)=XG(J2,J1) 
310:             ENDDO 
311:          ENDDO  
312:       ELSE 
313:          DO J1=1,N 
314:             MJ1=MOD(J1,NADDTARGET)+1 
315:             J3=3*(J1-1) 
316:             DO J2=J1+1,N 
317:                MJ2=MOD(J2,NADDTARGET)+1 
318:                J4=3*(J2-1) 
319:                R2T=(X(J3+1)-X(J4+1))**2+(X(J3+2)-X(J4+2))**2+(X(J3+3)-X(J4+3))**2 
320:                R2T=1.0D0/R2T 
321:                R6=R2T**3 
322:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
323:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*R6 
324:                ELSE 
325:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0) 
326:                ENDIF 
327:             ENDDO 
328:          ENDDO 
329:  
330:       ENDIF 
331:       ENERGY=4.0D0*ENERGY 
332:  
333:       IF (.NOT.GTEST) RETURN 
334:       DO J1=1,N 
335:          J3=3*J1 
336:          DUMMYX=0.0D0 
337:          DUMMYY=0.0D0 
338:          DUMMYZ=0.0D0 
339:          DO J4=1,N 
340:             J2=3*J4 
341:             XMUL2=XG(J4,J1) 
342:             DUMMYX=DUMMYX+XMUL2*(X(J3-2)-X(J2-2)) 
343:             DUMMYY=DUMMYY+XMUL2*(X(J3-1)-X(J2-1)) 
344:             DUMMYZ=DUMMYZ+XMUL2*(X(J3)  -X(J2)) 
345:          ENDDO 
346:          V(J3-2)=DUMMYX 
347:          V(J3-1)=DUMMYY 
348:          V(J3)=DUMMYZ 
349:       ENDDO 
350:  
351:       IF (.NOT.STEST) RETURN 
352:       CALL LJADDS(G,F,R2,R14,R8,X,N) 
353:  
354:       RETURN 
355:       END 
356:  
357: C***************************************************************************** 
358:  
359:       SUBROUTINE LJADDS2(G,F,R2,R14,R8,X,N) 
360:       USE MODHESS 
361:       USE COMMONS, ONLY : LJADDEPS, NADDTARGET 
362:       IMPLICIT NONE 
363:       INTEGER N, J1, J2, J3, J4, J5, J6, MJ1, MJ2 
364:       DOUBLE PRECISION G(N,N), R14(N,N), R8(N,N), 
365:      1                 R2(N,N), F(N,N),  
366:      2                 X(3*N),DUMMY 
367:  
368: C 
369: C  Calculate the g tensor. 
370: C 
371:       DO J1=1,N 
372:          MJ1=MOD(J1,NADDTARGET)+1 
373:          G(J1,J1)=0.0D0 
374:          DO J2=J1+1,N 
375:             MJ2=MOD(J2,NADDTARGET)+1 
376:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
377:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R14(J2,J1) 
378:             ELSE 
379:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R14(J2,J1)-R8(J2,J1)) 
380:             ENDIF 
381:             G(J1,J2)=G(J2,J1) 
382:          ENDDO 
383:       ENDDO 
384:  
385:       DO J1=1,N 
386:          MJ1=MOD(J1,NADDTARGET)+1 
387:          F(J1,J1)=0.0D0 
388:          DO J2=J1+1,N  
389:             MJ2=MOD(J2,NADDTARGET)+1 
390:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only 
391:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1) 
392:             ELSE 
393:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1)-LJADDEPS(MJ2,MJ1)*192.0D0*R8(J2,J1) 
394:             ENDIF 
395:             F(J1,J2)=F(J2,J1) 
396:          ENDDO 
397:       ENDDO 
398: C 
399: C  Now do the hessian. First are the entirely diagonal terms. 
400: C 
401:       DO J1=1,N 
402:          DO J2=1,3 
403:             J3=3*(J1-1)+J2 
404:             DUMMY=0.0D0 
405:             DO J4=1,N 
406:                DUMMY=DUMMY+F(J4,J1)*R2(J4,J1)* 
407:      1                 (X(J3)-X(3*(J4-1)+J2))**2 + G(J4,J1)    
408:             ENDDO 
409:             HESS(J3,J3)=DUMMY 
410:          ENDDO 
411:       ENDDO 
412: C 
413: C  Next are the terms where x_i and x_j are on the same atom 
414: C  but are different, e.g. y and z. 
415: C 
416:       DO J1=1,N 
417:          DO J2=1,3 
418:             J3=3*(J1-1)+J2 
419:             DO J4=J2+1,3 
420:                DUMMY=0.0D0 
421:                DO J5=1,N 
422:                   DUMMY=DUMMY + F(J5,J1)*R2(J5,J1)*  
423:      1           (X(J3)-X(3*(J5-1)+J2))*(X(3*(J1-1)+J4)-X(3*(J5-1)+J4))  
424:                ENDDO 
425:                HESS(3*(J1-1)+J4,J3)=DUMMY 
426:             ENDDO 
427:          ENDDO 
428:       ENDDO 
429: C 
430: C  Case III, different atoms, same cartesian coordinate. 
431: C 
432:       DO J1=1,N 
433:          DO J2=1,3 
434:             J3=3*(J1-1)+J2 
435:             DO J4=J1+1,N 
436:                HESS(3*(J4-1)+J2,J3)=-F(J4,J1)*R2(J4,J1)* 
437:      1                           (X(J3)-X(3*(J4-1)+J2))**2-G(J4,J1)  
438:             ENDDO 
439:          ENDDO 
440:       ENDDO 
441: C 
442: C  Case IV: different atoms and different cartesian coordinates. 
443: C 
444:       DO J1=1,N 
445:          DO J2=1,3 
446:             J3=3*(J1-1)+J2 
447:             DO J4=J1+1,N 
448:                DO J5=1,J2-1 
449:                   J6=3*(J4-1)+J5 
450:                   HESS(J6,J3)=-F(J4,J1)*R2(J4,J1) 
451:      1                    *(X(J3)-X(3*(J4-1)+J2)) 
452:      2                    *(X(3*(J1-1)+J5)-X(J6)) 
453:                   HESS(3*(J4-1)+J2,3*(J1-1)+J5)=HESS(J6,J3) 
454:                ENDDO 
455:             ENDDO 
456:          ENDDO 
457:       ENDDO 
458: C 
459: C  Symmetrise Hessian 
460: C 
461:       DO J1=1,3*N 
462:          DO J2=J1+1,3*N 
463:             HESS(J1,J2)=HESS(J2,J1) 
464:          ENDDO 
465:       ENDDO 
466:       RETURN 
467:       END 
468:  
469:  


r31389/potential.f90 2017-01-21 10:35:33.586250000 +0000 r31388/potential.f90 2017-01-21 10:35:34.214250000 +0000
168:       IF (FTEST) THEN168:       IF (FTEST) THEN
169:          RETURN169:          RETURN
170:       END IF170:       END IF
171:    ELSE IF (FRAUSIT) THEN171:    ELSE IF (FRAUSIT) THEN
172:       CALL RAD(X, GRAD, EREAL, GRADT)172:       CALL RAD(X, GRAD, EREAL, GRADT)
173:       CALL FRAUSI(NATOMS, X, GRAD, EREAL, GRADT, ANGST, NATOMS)173:       CALL FRAUSI(NATOMS, X, GRAD, EREAL, GRADT, ANGST, NATOMS)
174:       IF (FTEST) THEN174:       IF (FTEST) THEN
175:          RETURN175:          RETURN
176:       END IF176:       END IF
177:    ELSE IF (LJADDT) THEN177:    ELSE IF (LJADDT) THEN
178: 178:       CALL LJADD(NATOMS, X, GRAD, EREAL, GRADT, SECT)
179:       IF (RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN 
180:           XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS) 
181:           CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS) 
182:       ENDIF 
183:  
184:       IF (LJADD2T) THEN 
185:          CALL LJADD2(NATOMS, X, GRAD, EREAL, GRADT, SECT) 
186:       ELSE 
187:          CALL LJADD(NATOMS, X, GRAD, EREAL, GRADT, SECT) 
188:       ENDIF 
189:  
190:       IF ( RIGIDINIT .AND. (.NOT. ATOMRIGIDCOORDT)) THEN 
191:  
192:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
193:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS) 
194:  
195:          IF(GRADT) THEN 
196:             CALL TRANSFORMGRAD(GRAD, XRIGIDCOORDS, XRIGIDGRAD) 
197:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
198:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS) 
199:          ENDIF 
200:  
201:          IF(SECT) THEN 
202:             CALL TRANSFORMHESSIAN(HESS, GRAD, XRIGIDCOORDS, XRIGIDHESS, RBAANORMALMODET) 
203:             HESS(DEGFREEDOMS+1:3*NATOMS,:)=0.0D0 
204:             HESS(:, DEGFREEDOMS+1:3*NATOMS)=0.0D0 
205:             HESS(1:DEGFREEDOMS, 1:DEGFREEDOMS)=XRIGIDHESS(1:DEGFREEDOMS, 1:DEGFREEDOMS) 
206:          ENDIF 
207:       ENDIF 
208:  
209: !179: !
210: !  DIM Ne^+, Ne*, Ar^+, Ar*180: !  DIM Ne^+, Ne*, Ar^+, Ar*
211: !181: !
212:    ELSE IF ((NEON .OR. ARGON) .AND. (PLUS .OR. STAR)) THEN182:    ELSE IF ((NEON .OR. ARGON) .AND. (PLUS .OR. STAR)) THEN
213:       CALL RAD(X, GRAD, EREAL, GRADT)183:       CALL RAD(X, GRAD, EREAL, GRADT)
214: !      CALL RGNI(NATOMS, X, GRAD, EREAL, GRADT, h0, h1, ee, ev, w, NATOMS)184: !      CALL RGNI(NATOMS, X, GRAD, EREAL, GRADT, h0, h1, ee, ev, w, NATOMS)
215:       CALL RGNI(NATOMS, X, GRAD, EREAL, GRADT)185:       CALL RGNI(NATOMS, X, GRAD, EREAL, GRADT)
216: !186: !
217: !  DIM Ar^{2+}187: !  DIM Ar^{2+}
218: !188: !


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0