Revision dff3b5c68673d4d4c4c9f54c8ba110b8416098f2 authored by nitadori on 31 May 2020, 13:04:00 UTC, committed by GitHub on 31 May 2020, 13:04:00 UTC
Sorry, in the last commit I modified wrong part. In case !(JCL>N), SEMIX can be NaN.
1 parent 0b3202d
chpot.f
SUBROUTINE CHPOT(DP)
*
*
* Potential energy correction due to chain.
* -----------------------------------------
*
INCLUDE 'common6.h'
COMMON/CHAINC/ XC(3,NCMAX),UC(3,NCMAX),BODYC(NCMAX),ICH,
& LISTC(LMAX)
*
*
* Obtain current global coordinates of chain components.
CALL XCPRED(0)
*
* Consider contributions from all active perturbers.
DP = 0.0D0
JDUM = 0
NNBC = LISTC(1) + 1
*
DO 10 L = 2,NNBC
* Subtract potential energy due to chain c.m.
J = LISTC(L)
* Replace any regularized c.m. body by individual components.
IF (J.GT.N) THEN
JDUM = 2*(J - N) - 1
J = JDUM
END IF
2 A1 = X(1,ICH) - X(1,J)
A2 = X(2,ICH) - X(2,J)
A3 = X(3,ICH) - X(3,J)
RIJ2 = A1*A1 + A2*A2 + A3*A3
DP = DP - BODY(ICH)*BODY(J)/SQRT(RIJ2)
*
* Add individual interactions to obtain differential correction.
DO 5 K = 1,NCH
A1 = XC(1,K) - X(1,J)
A2 = XC(2,K) - X(2,J)
A3 = XC(3,K) - X(3,J)
RIJ2 = A1*A1 + A2*A2 + A3*A3
DP = DP + BODYC(K)*BODY(J)/SQRT(RIJ2)
5 CONTINUE
*
* Check for possible second KS component and restore dummy index.
IF (J.EQ.JDUM) THEN
J = JDUM + 1
JDUM = 0
GO TO 2
END IF
10 CONTINUE
*
RETURN
*
END
Computing file changes ...