X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc%2Fwham_calc1.F;h=ce5de230423b57da17f79d3b5a10a1a1747d1399;hb=aa037b39c46474f358b050551261d7c5ac3d2357;hp=676d4f4293e6bf2feef1cb18cfb7ffa06a254063;hpb=56b6e687a929b4935fe39b4abd8beb4e05265199;p=unres.git diff --git a/source/wham/src/wham_calc1.F b/source/wham/src/wham_calc1.F index 676d4f4..ce5de23 100644 --- a/source/wham/src/wham_calc1.F +++ b/source/wham/src/wham_calc1.F @@ -84,7 +84,9 @@ c parameter (MaxHdim=200000) & 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 @@ -94,7 +96,6 @@ c parameter (MaxHdim=200000) 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", @@ -220,8 +221,8 @@ c parameter (MaxHdim=200000) 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 @@ -306,10 +307,16 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft 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 @@ -321,7 +328,8 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft & +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) @@ -331,7 +339,8 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft & +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), @@ -347,6 +356,24 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft 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 @@ -361,6 +388,8 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft & etot,v(i,kk,ib,iparm) #endif enddo ! kk + + endif enddo ! ib enddo ! iparm enddo ! i @@ -580,7 +609,7 @@ 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) + & (enetb(k,i,iparm),k=1,22) #endif call restore_parm(iparm) #ifdef DEBUG @@ -665,10 +694,19 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft 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 @@ -680,7 +718,8 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft & +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) @@ -690,14 +729,19 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft & +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 -c write (iout,*) "i",i," ib",ib, -c & " temp",1.0d0/(1.987d-3*beta_h(ib,iparm))," etot",etot, -c & " entfac",entfac(i) +#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) etot=etot-entfac(i)/beta_h(ib,iparm) +#endif 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 @@ -910,6 +954,13 @@ c write (iout,*) "me1",me1," scount",scount(me1) 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) do k=0,nGridT betaT=startGridT+k*delta_T temper=betaT @@ -1018,7 +1069,7 @@ c write (iout,*) ib," PotEmin",potEmin & +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 eprim=ftprim(6)*evdw_t+ftprim(1)*welec*ees & +ftprim(1)*wtor*etors+ & ftprim(3)*wcorr*ecorr+ftprim(4)*wcorr5*ecorr5+ @@ -1041,7 +1092,7 @@ c write (iout,*) ib," PotEmin",potEmin & +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 eprim=ftprim(6)*evdw_t+ftprim(1)*welec*(ees+evdw1) & +ftprim(1)*wtor*etors+ & ftprim(3)*wcorr*ecorr+ftprim(4)*wcorr5*ecorr5+ @@ -1058,6 +1109,11 @@ c write (iout,*) ib," PotEmin",potEmin #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,