& eplus,eminus,logfac,tanhT,tt
double precision etot,evdw,evdw_t,evdw2,ees,evdw1,ebe,etors,
& escloc,ehpb,ecorr,ecorr5,ecorr6,eello_turn4,eello_turn3,
- & eturn6,eel_loc,edihcnstr,etors_d,estr,evdw2_14,esccor
+ & eturn6,eel_loc,edihcnstr,etors_d,estr,evdw2_14,esccor,
+ & ehomology_constr,edfadis,edfator,edfanei,edfabet
+
integer ind_point(maxpoint),upindE,indE
character*16 plik
character*128 nazwa
integer ilen
external ilen
-
write(licz2,'(bz,i2.2)') islice
nbin1 = 1.0d0/delta
write (iout,'(//80(1h-)/"Solving WHAM equations for slice",
c write (9,'(3i5,f10.5)') i,(iparm,potE(i,iparm),iparm=1,nParmSet)
do iparm=1,nParmSet
#ifdef DEBUG
- write (iout,'(2i5,21f8.2)') i,iparm,
- & (enetb(k,i,iparm),k=1,21)
+ write (iout,'(2i5,22f8.2)') i,iparm,
+ & (enetb(k,i,iparm),k=1,22)
#endif
call restore_parm(iparm)
#ifdef DEBUG
estr=enetb(18,i,iparm)
esccor=enetb(19,i,iparm)
edihcnstr=enetb(20,i,iparm)
+ ehomology_constr=enetb(22,i,iparm)
+ edfadis=enetb(23,i,iparm)
+ edfator=enetb(24,i,iparm)
+ edfanei=enetb(25,i,iparm)
+ edfabet=enetb(26,i,iparm)
#ifdef DEBUG
- write (iout,'(3i5,6f5.2,14f12.3)') i,ib,iparm,(ft(l),l=1,6),
+ write (iout,'(3i5,6f5.2,16f12.3)') i,ib,iparm,(ft(l),l=1,6),
& evdw+evdw_t,evdw2,ees,evdw1,ecorr,eel_loc,estr,ebe,escloc,
- & etors,etors_d,eello_turn3,eello_turn4,esccor
+ & etors,etors_d,eello_turn3,eello_turn4,esccor,
+ & ehomology_constr,edfadis,edfator,edfanei,edfabet
#endif
#ifdef SPLITELE
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+ft(1)*welec*ees
& +wvdwpp*evdw1
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc
& +edihcnstr+ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
#else
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2
& +ft(1)*welec*(ees+evdw1)
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc+edihcnstr
& +ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
#endif
#ifdef DEBUG
write (iout,*) i,iparm,1.0d0/(beta_h(ib,iparm)*1.987D-3),
call enerprint(energia(0),fT)
endif
#endif
+#ifdef DEBUG
+ write (iout,*) "homol_nset",homol_nset,nR(ib,iparm)
+#endif
+ if (homol_nset.gt.1) then
+
+ do kk=1,nR(ib,iparm)
+ Econstr=waga_homology(kk)*ehomology_constr
+ v(i,kk,ib,iparm)=
+ & -beta_h(ib,iparm)*(etot+Econstr)
+#ifdef DEBUG
+ write (iout,'(4i5,4e15.5)') i,kk,ib,iparm,
+ & etot,Econstr,v(i,kk,ib,iparm)
+#endif
+ enddo ! kk
+
+ else
+
+ etot=etot+ehomology_constr
do kk=1,nR(ib,iparm)
Econstr=0.0d0
do j=1,nQ
& etot,v(i,kk,ib,iparm)
#endif
enddo ! kk
+
+ endif
enddo ! ib
enddo ! iparm
enddo ! i
do i=1,nR(ib,iparm)
fi(i,ib,iparm)=fi(i,ib,iparm)-avefi
enddo
- write (iout,'(8f10.5)') (fi(i,ib,iparm),i=1,nR(ib,iparm))
- write (iout,'(8f10.5)') (f(i,ib,iparm),i=1,nR(ib,iparm))
+ write (iout,'(6f15.5)') (fi(i,ib,iparm),i=1,nR(ib,iparm))
+ write (iout,'(6f15.5)') (f(i,ib,iparm),i=1,nR(ib,iparm))
enddo
enddo
do iparm=1,nParmSet
#ifdef DEBUG
write (iout,'(2i5,21f8.2)') i,iparm,
- & (enetb(k,i,iparm),k=1,21)
+ & (enetb(k,i,iparm),k=1,22)
#endif
call restore_parm(iparm)
#ifdef DEBUG
estr=enetb(18,i,iparm)
esccor=enetb(19,i,iparm)
edihcnstr=enetb(20,i,iparm)
+ ehomology_constr=enetb(22,i,iparm)
+ if (homol_nset.gt.1)
+ & ehomology_constr=waga_homology(ihset)*ehomology_constr
+ edfadis=enetb(23,i,iparm)
+ edfator=enetb(24,i,iparm)
+ edfanei=enetb(25,i,iparm)
+ edfabet=enetb(26,i,iparm)
#ifdef DEBUG
- write (iout,'(3i5,6f5.2,14f12.3)') i,ib,iparm,(ft(l),l=1,6),
+ write (iout,'(3i5,6f5.2,17f12.3)') i,ib,iparm,(ft(l),l=1,6),
& evdw+evdw_t,evdw2,ees,evdw1,ecorr,eel_loc,estr,ebe,escloc,
- & etors,etors_d,eello_turn3,eello_turn4,esccor,edihcnstr
+ & etors,etors_d,eello_turn3,eello_turn4,esccor,edihcnstr,
+ & ehomology_constr+wdfa_dist*edfadis+wdfa_tor*edfator+
+ & wdfa_nei*edfanei+wdfa_beta*edfabet
#endif
#ifdef SPLITELE
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+ft(1)*welec*ees
& +wvdwpp*evdw1
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc
& +edihcnstr+ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+ehomology_constr+wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
#else
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2
& +ft(1)*welec*(ees+evdw1)
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc+edihcnstr
& +ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+ehomology_constr+wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
+#endif
+#ifdef DEBUG
+ write (iout,*) "i",i," ib",ib,
+ & " temp",1.0d0/(1.987d-3*beta_h(ib,iparm))," etot",etot,
+ & " entfac",entfac(i)," ecorr",etot-entfac(i)/beta_h(ib,iparm)
#endif
-c write (iout,*) "i",i," ib",ib,
-c & " temp",1.0d0/(1.987d-3*beta_h(ib,iparm))," etot",etot,
-c & " entfac",entfac(i)
etot=etot-entfac(i)/beta_h(ib,iparm)
if(etot.lt.potEmin_all(ib,iparm)) potEmin_all(ib,iparm)=etot
-c write (iout,*) "efree",etot," potEmin",potEmin_all(ib,iparm)
+#ifdef DEBUG
+ write (iout,*) "efree",etot," potEmin",potEmin_all(ib,iparm)
+#endif
enddo ! ib
enddo ! iparm
enddo ! i
do i=1,nParmSet
write (iout,*) "Parameter set",i
do j=1,nT_h(i)
- write (iout,*) j,PotEmin_all(j,i)
+ write (iout,*) j,1.0d0/(1.987d-3*beta_h(j,i)),
+ & PotEmin_all(j,i)
enddo
enddo
write (iout,*) "potEmin_min",potEmin_min
estr=enetb(18,t,iparm)
esccor=enetb(19,t,iparm)
edihcnstr=enetb(20,t,iparm)
+ ehomology_constr=enetb(22,t,iparm)
+ if (homol_nset.gt.1)
+ & ehomology_constr=waga_homology(ihset)*ehomology_constr
+ edfadis=enetb(23,t,iparm)
+ edfator=enetb(24,t,iparm)
+ edfanei=enetb(25,t,iparm)
+ edfabet=enetb(26,t,iparm)
+c do k=0,nGridT
do k=0,nGridT
betaT=startGridT+k*delta_T
temper=betaT
c write (iout,*) "ftbis",ftbis
betaT=1.0d0/(1.987D-3*betaT)
if (betaT.ge.beta_h(1,iparm)) then
- potEmin=potEmin_all(1,iparm)
-c write(iout,*) "first",temper,potEmin
- else if (betaT.lt.beta_h(nT_h(iparm),iparm)) then
- potEmin=potEmin_all(nT_h(iparm),iparm)
-c write (iout,*) "last",temper,potEmin
+ potEmin=potEmin_all(1,iparm)+
+ & (potEmin_all(1,iparm)-potEmin_all(2,iparm))/
+ & (1.0/beta_h(1,iparm)-1.0/beta_h(2,iparm))*
+ & (1.0/betaT-1.0/beta_h(1,iparm))
+#ifdef DEBUG
+ write(iout,*) "first",temper,potEmin
+#endif
+ else if (betaT.le.beta_h(nT_h(iparm),iparm)) then
+ potEmin=potEmin_all(nT_h(iparm),iparm)+
+ &(potEmin_all(nT_h(iparm),iparm)-potEmin_all(nT_h(iparm)-1,iparm))/
+ &(1.0/beta_h(nT_h(iparm),iparm)-1.0/beta_h(nT_h(iparm)-1,iparm))*
+ &(1.0/betaT-1.0/beta_h(nt_h(iparm),iparm))
+#ifdef DEBUG
+ write (iout,*) "last",temper,potEmin
+#endif
else
do l=1,nT_h(iparm)-1
if (betaT.le.beta_h(l,iparm) .and.
& betaT.gt.beta_h(l+1,iparm)) then
potEmin=potEmin_all(l,iparm)
-c write (iout,*) "l",l,
-c & betaT,1.0d0/(1.987D-3*beta_h(l,iparm)),
-c & 1.0d0/(1.987D-3*beta_h(l+1,iparm)),temper,potEmin
+#ifdef DEBUG
+ write (iout,*) "l",l,
+ & betaT,1.0d0/(1.987D-3*beta_h(l,iparm)),
+ & 1.0d0/(1.987D-3*beta_h(l+1,iparm)),temper,potEmin
+#endif
exit
endif
enddo
endif
-c write (iout,*) ib," PotEmin",potEmin
+#ifdef DEBUG
+ write (iout,*) "k",k," potEmin",potEmin
+#endif
#ifdef SPLITELE
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+ft(1)*welec*ees
& +wvdwpp*evdw1
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc
& +edihcnstr+ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+ehomology_constr
+ & +wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
eprim=ftprim(6)*evdw_t+ftprim(1)*welec*ees
& +ftprim(1)*wtor*etors+
& ftprim(3)*wcorr*ecorr+ftprim(4)*wcorr5*ecorr5+
etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2
& +ft(1)*welec*(ees+evdw1)
& +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
- & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
+ & +wstrain*ehpb+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
& +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
& +ft(2)*wturn3*eello_turn3
& +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc+edihcnstr
& +ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
- & +wbond*estr
+ & +wbond*estr+ehomology_constr
+ & +wdfa_dist*edfadis
+ & +wdfa_tor*edfator+wdfa_nei*edfanei+wdfa_beta*edfabet
eprim=ftprim(6)*evdw_t+ftprim(1)*welec*(ees+evdw1)
& +ftprim(1)*wtor*etors+
& ftprim(3)*wcorr*ecorr+ftprim(4)*wcorr5*ecorr5+
#endif
weight=dexp(-betaT*(etot-potEmin)+entfac(t))
#ifdef DEBUG
+ write (iout,'(3i5,6f5.2,17f12.3)') i,ib,iparm,(ft(l),l=1,6),
+ & evdw+evdw_t,evdw2,ees,evdw1,ecorr,eel_loc,estr,ebe,escloc,
+ & etors,etors_d,eello_turn3,eello_turn4,esccor,edihcnstr,
+ & ehomology_constr+wdfa_dist*edfadis+wdfa_tor*edfator+
+ & wdfa_nei*edfanei+wdfa_beta*edfabet
write (iout,*) "iparm",iparm," t",t," temper",temper,
& " etot",etot," entfac",entfac(t),
& " efree",etot-entfac(t)/betaT," potEmin",potEmin,
else
write (iout,'(a,i3)') "Parameter set",iparm
endif
+c do i=0,NGridT
do i=0,NGridT
betaT=1.0d0/(1.987D-3*(startGridT+i*delta_T))
if (betaT.ge.beta_h(1,iparm)) then
endif
enddo
endif
+
+c write (iout,*) "i",i," potEmin",potEmin
+
sumE(i,iparm)=sumE(i,iparm)/sumW(i,iparm)
sumEbis(i,iparm)=(startGridT+i*delta_T)*sumEbis(i,iparm)/
& sumW(i,iparm)
projects / unres.git / blobdiff