X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc-M%2Fenergy_p_new.F;h=fb18913435cc4daf46f956bdf7c2a19b66cf1132;hb=48f04f24e913a3e10867d2038b30efcd48a60a9f;hp=2b993941a9ab62915800611584d4f7da4408d0e2;hpb=d305127baa0f7e25ecb5422e6e0e34aa15db4776;p=unres.git diff --git a/source/wham/src-M/energy_p_new.F b/source/wham/src-M/energy_p_new.F index 2b99394..fb18913 100644 --- a/source/wham/src-M/energy_p_new.F +++ b/source/wham/src-M/energy_p_new.F @@ -224,6 +224,11 @@ c detecting NaNQ #ifdef MPL c endif #endif +#define DEBUG +#ifdef DEBUG + call enerprint(energia,fact) +#endif +#undef DEBUG if (calc_grad) then C C Sum up the components of the Cartesian gradient. @@ -265,12 +270,29 @@ C & wturn6*fact(5)*gcorr6_turn(j,i)+ & wsccor*fact(2)*gsccorc(j,i) & +wliptran*gliptranc(j,i) + & +welec*gshieldc(j,i) + & +welec*gshieldc_loc(j,i) + & +wcorr*gshieldc_ec(j,i) + & +wcorr*gshieldc_loc_ec(j,i) + & +wturn3*gshieldc_t3(j,i) + & +wturn3*gshieldc_loc_t3(j,i) + & +wturn4*gshieldc_t4(j,i) + & +wturn4*gshieldc_loc_t4(j,i) + & +wel_loc*gshieldc_ll(j,i) + & +wel_loc*gshieldc_loc_ll(j,i) + gradx(j,i,icg)=fact(1)*wsc*gvdwx(j,i) & +fact(1)*wscp*gradx_scp(j,i)+ & wbond*gradbx(j,i)+ & wstrain*ghpbx(j,i)+wcorr*gradxorr(j,i)+ & wsccor*fact(2)*gsccorx(j,i) & +wliptran*gliptranx(j,i) + & +welec*gshieldx(j,i) + & +wcorr*gshieldx_ec(j,i) + & +wturn3*gshieldx_t3(j,i) + & +wturn4*gshieldx_t4(j,i) + & +wel_loc*gshieldx_ll(j,i) + endif enddo @@ -309,12 +331,29 @@ C & wturn6*fact(5)*gcorr6_turn(j,i)+ & wsccor*fact(2)*gsccorc(j,i) & +wliptran*gliptranc(j,i) + & +welec*gshieldc(j,i) + & +welec*gshieldc_loc(j,i) + & +wcorr*gshieldc_ec(j,i) + & +wcorr*gshieldc_loc_ec(j,i) + & +wturn3*gshieldc_t3(j,i) + & +wturn3*gshieldc_loc_t3(j,i) + & +wturn4*gshieldc_t4(j,i) + & +wturn4*gshieldc_loc_t4(j,i) + & +wel_loc*gshieldc_ll(j,i) + & +wel_loc*gshieldc_loc_ll(j,i) + gradx(j,i,icg)=fact(1)*wsc*gvdwx(j,i)+ & fact(1)*wscp*gradx_scp(j,i)+ & wbond*gradbx(j,i)+ & wstrain*ghpbx(j,i)+wcorr*gradxorr(j,i)+ & wsccor*fact(1)*gsccorx(j,i) & +wliptran*gliptranx(j,i) + & +welec*gshieldx(j,i) + & +wcorr*gshieldx_ec(j,i) + & +wturn3*gshieldx_t3(j,i) + & +wturn4*gshieldx_t4(j,i) + & +wel_loc*gshieldx_ll(j,i) + endif enddo #endif @@ -1034,7 +1073,12 @@ C lipbufthick is thickenes of lipid buffore & +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 bb=bb_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & +bb_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 -C write(iout,*) "tu,", i,j,aa_lip(itypi,itypj),bb_lip(itypi,itypj) +C if (aa.ne.aa_aq(itypi,itypj)) then + +C write(iout,*) "tu,", i,j,aa_aq(itypi,itypj)-aa, +C & bb_aq(itypi,itypj)-bb, +C & sslipi,sslipj +C endif C write(iout,*),aa,aa_lip(itypi,itypj),aa_aq(itypi,itypj) C checking the distance @@ -1118,6 +1162,7 @@ c & aux*e2/eps(itypi,itypj) c if (lprn) then sigm=dabs(aa/bb)**(1.0D0/6.0D0) epsi=bb**2/aa +C#define DEBUG #ifdef DEBUG write (iout,'(2(a3,i3,2x),17(0pf7.3))') & restyp(itypi),i,restyp(itypj),j, @@ -1127,6 +1172,7 @@ c if (lprn) then & evdwij write (iout,*) "partial sum", evdw, evdw_t #endif +C#undef DEBUG c endif if (calc_grad) then C Calculate gradient components. @@ -2227,6 +2273,12 @@ c write (iout,*) "i",i,iteli," j",j,itelj," eesij",eesij C 12/26/95 - for the evaluation of multi-body H-bonding interactions ees0ij=4.0D0+fac*fac-3.0D0*(cosb*cosb+cosg*cosg) if (shield_mode.gt.0) then +C#define DEBUG +#ifdef DEBUG + write(iout,*) "ees_compon",i,j,el1,el2, + & fac_shield(i),fac_shield(j) +#endif +C#undef DEBUG C fac_shield(i)=0.4 C fac_shield(j)=0.6 el1=el1*fac_shield(i)**2*fac_shield(j)**2 @@ -2962,6 +3014,8 @@ C Derivatives due to the contact function & *fac_shield(i)*fac_shield(j) gacontp_hb3(k,num_conti,i)=gggp(k) + & *fac_shield(i)*fac_shield(j) + gacontm_hb1(k,num_conti,i)=ghalfm & +(ecosam*(dc_norm(k,j)-cosa*dc_norm(k,i)) & + ecosbm*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1) @@ -3040,7 +3094,7 @@ C end of changes suggested by Ana C & .or. itype(i+5).eq.ntyp1 C & .or. itype(i).eq.ntyp1 C & .or. itype(i-1).eq.ntyp1 - & ) goto 178 + & ) goto 179 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C @@ -3121,6 +3175,7 @@ C Derivatives in gamma(i) call transpose2(auxmat2(1,1),pizda(1,1)) call matmat2(a_temp(1,1),pizda(1,1),pizda(1,1)) gel_loc_turn3(i)=gel_loc_turn3(i)+0.5d0*(pizda(1,1)+pizda(2,2)) + & *fac_shield(i)*fac_shield(j) C Derivatives in gamma(i+1) call matmat2(EUg(1,1,i+1),EUgder(1,1,i+2),auxmat2(1,1)) call transpose2(auxmat2(1,1),pizda(1,1)) @@ -3169,6 +3224,7 @@ C Cartesian derivatives enddo endif + 179 continue else if (j.eq.i+3 .and. itype(i+2).ne.ntyp1) then if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1 C changes suggested by Ana to avoid out of bounds