*Deck Prop1E line.1 Subroutine Prop1E(V,MDV) line.2 C The options which control which properties are to be calculated line.77 C are described in the comments in L601. line.78 C line.79 C%SN Kosugi's surface (1993) was incorpolated by S. Nakagawa < insert 6 lines C%SN Electro Static Potentials are written to REST file < comment C%SN 2001 Jan 26 for G98 < C%SN 2005 Feb 27 for G03 < C%SN 2015 Nov 28 for G09 < C%SN < #include "commonmol.inc" line.80 c.............................................................................................. Real*8 V(*),OrigGd(3,MaxGrd),StepGd(3,3,MaxGrd),ScaDFX(4), line.94 $ D1P(3),D2P(6),Dipole(3),CelVec(3,3),XX(1) > change to this line C%SN $ D1P(3),D2P(6),Dipole(3),XX(1) line.95 Integer IDerGd(2,MaxGrd),NPtSid(3,MaxGrd),KTape(MaxGrd), line.96 c............................................................................................. Call AddBMV(IRwMul,1,NAtoms,V(IESPCh),V(jAvail),NAvDB) line.625 C%SN< < insert 6 lines if(IGrTyp.ge.1.and.IGrTyp.le.9) then < CALL WriESP(IOut,NAtoms,NFit,V(jElect),V(jFit),IAn,C, < $ V(IESPCh),IDens) < endif < C%SN> < If(WrtFit) Call WtESPF(IOut,IPrint,ICharg,Multip,NAtoms,IAn, line.626 $ C,V(IESPCh),NFit,V(jFit),V(jElect+NAtoms),V(jAvail), line.627 $ NAvail) line.628 c---------------------------------------------------------------------------------------------- *Deck ChgGrd line.961 Subroutine ChgGrd(In,IOut,IPrint,IDump,IGrTyp,IVDTyp,IRedVD, line.962 else if(IGrTyp.eq.4) then line.1032 Call GPGrid(IOut,IPrint,NAtUse,IAnL,CL,Rad,MxCent,NDens, line.1033 $ LayerI,SIncIn,NFit,Cent,V,MDV) line.1034 C%SN< < insert 22 lines else if(IGrTyp.eq.6) then < write(IOut,*) ' Kosugi vdw 2.0 surface' < ESFACT=2.0d0 < Call SurGS3(IOut,IPrint,NAtoms,C,Rad,ESFACT,NFit, < $ Cent) < else if(IGrTyp.eq.7) then < write(IOut,*) ' Kosugi vdw 1.7 surface' < ESFACT=1.7d0 < Call SurGS3(IOut,IPrint,NAtoms,C,Rad,ESFACT,NFit, < $ Cent) < else if(IGrTyp.eq.8) then < write(IOut,*) ' Kosugi vdw 1.8 surface' < ESFACT=1.8d0 < Call SurGS3(IOut,IPrint,NAtoms,C,Rad,ESFACT,NFit, < $ Cent) < else if(IGrTyp.eq.9) then < write(IOut,*) ' Kosugi vdw 1.9 surface' < ESFACT=1.9d0 < Call SurGS3(IOut,IPrint,NAtoms,C,Rad,ESFACT,NFit, < $ Cent) < C%SN> < else line.1035 Write(IOut,1000) IGrTyp line.1036 c---------------------------------------------------------------------------------------------- *Deck ESPFit line.1257 Subroutine ESPFit(IOut,IPrint,ToAng,ToE,FitAll,FitDip,UseDip, line.1258 C TotChg = total charge line.1275 C%SN efg(ESP) = TOTAL ESP AT THE POINTS < insert 3 lines C%SN cent(POTPT) = ESP POINTS < comment C%SN C(CO) = COORDINATES < C EFG = ESP at the points line.1276 c---------------------------------------------------------------------------------------------- *Deck SurGS3 << Add Subroutine C%SN< << SurGS3 to the bottom Subroutine SurGS3(IOut,IPrint,NAtoms,CO,Rad,ESFACT,IESTAT,Cent) << line of l602.F IMPLICIT REAL*8(A-H,O-Z) Dimension CO(3,NAtoms),Rad(NAtoms),Cent(3,*) Dimension CR(3200),THETA(20) WRITE(IOut,*) ' ' IF(ESFACT.EQ.0.D0) ESFACT=1.8D0 RADIAN= DATAN(1.D0)*4.D0 TWORAD= RADIAN*2 HLFRAD= RADIAN*0.5D0 RAQ= RADIAN/10.D0 DO 100 I=1,11 100 THETA(I)= RAQ*(I-1)-HLFRAD DIV=TWORAD/4.D0 C%SN< do IAT=1,NAtoms CR(IAT)=Rad(IAT)*ESFACT/0.529177D0 enddo C%SN> IP=0 DO 1000 JAT=1,NAtoms IF(CR(JAT).LT.1.D0) GO TO 1000 C IJ=0 FCT=-1.D0 DO 200 I=1,11 FCT=-FCT RDCOS= CR(JAT)*DCOS(THETA(I))*FCT RDSCZ=CO(3,JAT)+CR(JAT)*DSIN(THETA(I)) DIA=TWORAD*DABS(RDCOS) RDIV=DIA/DIV IDIV=RDIV+1.0001D0 TWORAQ= TWORAD/IDIV DO 240 J=1,IDIV PHAI=TWORAQ*(J-1)-RADIAN IJ=IJ+1 Cent(1,IP+IJ)=CO(1,JAT)+RDCOS*DSIN(PHAI) Cent(2,IP+IJ)=CO(2,JAT)+RDCOS*DCOS(PHAI) 240 Cent(3,IP+IJ)=RDSCZ 200 CONTINUE C NP=IP DO 1100 I=IP+1,IP+IJ DO 1110 IAT=1,NAtoms IF(IAT.EQ.JAT.OR.CR(IAT).LT.1.D0) GO TO 1110 IF((Cent(1,I)-CO(1,IAT))**2+(Cent(2,I)-CO(2,IAT))**2 & +(Cent(3,I)-CO(3,IAT))**2-CR(IAT)**2.LT.-1.D-3) GO TO 1100 1110 CONTINUE NP=NP+1 Cent(1,NP)=Cent(1,I) Cent(2,NP)=Cent(2,I) Cent(3,NP)=Cent(3,I) 1100 CONTINUE NEIGH=IP+IJ-NP IESTAT=NP IF(IP.EQ.0) GO TO 1900 NP=IP C < C Compensate for nearby surface points C DO 1200 I=IP+1,IESTAT DO 1210 J=1,IP 1210 IF((Cent(1,I)-Cent(1,J))**2+(Cent(2,I)-Cent(2,J))**2 & +(Cent(3,I)-Cent(3,J))**2.LT.0.6D0*ESFACT) GO TO 1200 NP=NP+1 Cent(1,NP)=Cent(1,I) Cent(2,NP)=Cent(2,I) Cent(3,NP)=Cent(3,I) 1200 CONTINUE C > 1900 IP=NP WRITE(IOut,'(I10,f10.4,6X,A,I6,6X,A,I6)') * JAT,Rad(JAT),'NEIGHBOR CUT',NEIGH,'INTERSECT CUT',IESTAT-NP 1000 CONTINUE IESTAT=IP WRITE(IOut,'( 16X,A,I6 )') 'TOTAL POINTS',IESTAT WRITE(IOut,'(/16X,A,F10.3)') * 'HOW MANY TIMES THE VAN DER WAALS RADIUS',ESFACT do i=1,iestat do j=1,3 Cent(j,i)=Cent(j,i)*0.529177d0 enddo enddo RETURN END c--------------------------------------------------------------------------------------------------- *Deck WriESP << Add Subroutine Subroutine WriESP(IOut,NAtoms,NFit,ESP,Cent,IAn,C,Q,IDens) << WriESP to the bottom Implicit Real*8(A-H,O-Z) << line of Subroutine << SurGS3 character*8 rest Dimension IAn(NAtoms), C(3,NAtoms), Q(NAtoms), $ ESP(NFit,10), Cent(3,NFit) data au/0.529177d0/ IV=1 Charge=0.0d0 if(IDens.eq.0)then rest='Rest_scf' else if(IDens.eq.1)then rest='Rest_mp1' else if(IDens.eq.2)then rest='Rest_mp2' else if(IDens.eq.3)then rest='Rest_mp3' else if(IDens.eq.4)then rest='Rest_mp4' else if(IDens.eq.5)then rest='Rest_cio' else if(IDens.eq.6)then rest='Rest_cid' else if(IDens.eq.7)then rest='Rest_qci' else if(IDens.eq.8)then rest='Rest_m2s' else rest='Rest_xxx' endif Open(1,file=rest,Access='Sequential',Form='Formatted') Write(1,*) NAtoms,Charge Do i=1,NAtoms Write(1,'(I2,4(F16.10,I2))') & IAn(i),C(1,i),IV,C(2,i),IV,C(3,i),IV,Q(i),IV Enddo Write(1,*) NFit Do i=1,NFit Write(1,'(4f15.9)') & Cent(1,i),Cent(2,i),Cent(3,i),-ESP(i+NAtoms,1) Enddo CLOSE(Unit=1) RETURN END C%SN>