From: Adam Sieradzan Date: Tue, 22 Sep 2015 08:10:57 +0000 (+0200) Subject: working shielding X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?p=unres.git;a=commitdiff_plain;h=2095e643aafa2fb192e066e69e5ca3b1a59c2727 working shielding --- diff --git a/source/unres/src_MD-M/COMMON.DERIV b/source/unres/src_MD-M/COMMON.DERIV index 8082b0a..17e3184 100644 --- a/source/unres/src_MD-M/COMMON.DERIV +++ b/source/unres/src_MD-M/COMMON.DERIV @@ -2,7 +2,7 @@ & gvdwpp,gel_loc,gel_loc_long,gvdwc_scpp,gliptranc,gliptranx, & gradx_scp,gvdwc_scp,ghpbx,ghpbc,gloc,gloc_x,dtheta,dphi,dalpha, & domega,gscloc,gsclocx,gradcorr,gradcorr_long,gradcorr5_long, - & gradcorr6_long,gcorr6_turn_long,gvdwx + & gradcorr6_long,gcorr6_turn_long,gvdwx,gshieldx integer nfl,icg common /derivat/ dcdv(6,maxdim),dxdv(6,maxdim),dxds(6,maxres), & gradx(3,-1:maxres,2),gradc(3,-1:maxres,2),gvdwx(3,-1:maxres), @@ -10,6 +10,9 @@ & gvdwpp(3,-1:maxres),gvdwc_scpp(3,-1:maxres), & gliptranc(3,-1:maxres), & gliptranx(3,-1:maxres), + & gshieldx(3,-1:maxres), + & gshieldc(3,-1:maxres), + & gshieldc_loc(3,-1:maxres), & gradafm(3,-1:maxres), & gradx_scp(3,-1:maxres),gvdwc_scp(3,-1:maxres), & ghpbx(3,-1:maxres), diff --git a/source/unres/src_MD-M/COMMON.SHIELD b/source/unres/src_MD-M/COMMON.SHIELD index 7974844..65d7891 100644 --- a/source/unres/src_MD-M/COMMON.SHIELD +++ b/source/unres/src_MD-M/COMMON.SHIELD @@ -1,9 +1,12 @@ double precision VSolvSphere,VSolvSphere_div,long_r_sidechain, - & short_r_sidechain,fac_shield + & short_r_sidechain,fac_shield,grad_shield_side,grad_shield, + & buff_shield integer ishield_list,shield_list - common /shield/ VSolvSphere,VSolvSphere_div, + common /shield/ VSolvSphere,VSolvSphere_div,buff_shield, & long_r_sidechain(ntyp), & short_r_sidechain(ntyp),fac_shield(maxres), - & ishield_list(maxres),shield_list(maxres) + & grad_shield_side(3,15,-1:maxres),grad_shield(3,-1:maxres), + & grad_shield_loc(3,15,-1:maxres), + & ishield_list(maxres),shield_list(15,maxres) diff --git a/source/unres/src_MD-M/COMMON.TORSION b/source/unres/src_MD-M/COMMON.TORSION index 5515058..4cf97e3 100644 --- a/source/unres/src_MD-M/COMMON.TORSION +++ b/source/unres/src_MD-M/COMMON.TORSION @@ -27,14 +27,14 @@ C surface double precision b1,b2,cc,dd,ee,ctilde,dtilde,b2tilde,b1tilde, &bnew1,bnew2,eenew,gtb1,gtb2,eeold,gtee integer nloctyp - common/fourier/ b1(2,maxres),b2(2,maxres),b(13,0:maxtor), + common/fourier/ b1(2,0:maxres),b2(2,0:maxres),b(13,0:maxtor), & bnew1(3,2,-maxtor:maxtor),bnew2(3,2,-maxtor:maxtor), & cc(2,2,-maxtor:maxtor), & dd(2,2,-maxtor:maxtor),eeold(2,2,-maxtor:maxtor), & eenew(2,-maxtor:maxtor), & ee(2,2,maxres), & ctilde(2,2,-maxtor:maxtor), - & dtilde(2,2,-maxtor:maxtor),b1tilde(2,maxres), + & dtilde(2,2,-maxtor:maxtor),b1tilde(2,0:maxres), & b2tilde(2,maxres), & gtb1(2,maxres),gtb2(2,maxres),gtEE(2,2,maxres), & nloctyp diff --git a/source/unres/src_MD-M/energy_p_new_barrier.F b/source/unres/src_MD-M/energy_p_new_barrier.F index 109af35..5402971 100644 --- a/source/unres/src_MD-M/energy_p_new_barrier.F +++ b/source/unres/src_MD-M/energy_p_new_barrier.F @@ -141,7 +141,7 @@ C Introduction of shielding effect first for each peptide group C the shielding factor is set this factor is describing how each C peptide group is shielded by side-chains C the matrix - shield_fac(i) the i index describe the ith between i and i+1 - write (iout,*) "shield_mode",shield_mode +C write (iout,*) "shield_mode",shield_mode if (shield_mode.gt.0) then call set_shield_fac endif @@ -546,6 +546,7 @@ c enddo & wstrain*ghpbc(j,i) & +wliptran*gliptranc(j,i) & +gradafm(j,i) + & +welec*gshieldc(j,i) enddo enddo @@ -564,6 +565,7 @@ c enddo & wstrain*ghpbc(j,i) & +wliptran*gliptranc(j,i) & +gradafm(j,i) + & +welec*gshieldc(j,i) enddo enddo @@ -682,6 +684,13 @@ c enddo do i=-1,nct do j=1,3 #ifdef SPLITELE +C print *,gradbufc(1,13) +C print *,welec*gelc(1,13) +C print *,wel_loc*gel_loc(1,13) +C print *,0.5d0*(wscp*gvdwc_scpp(1,13)) +C print *,welec*gelc_long(1,13)+wvdwpp*gvdwpp(1,13) +C print *,wel_loc*gel_loc_long(1,13) +C print *,gradafm(1,13),"AFM" gradc(j,i,icg)=gradbufc(j,i)+welec*gelc(j,i)+ & wel_loc*gel_loc(j,i)+ & 0.5d0*(wscp*gvdwc_scpp(j,i)+ @@ -702,6 +711,10 @@ c enddo & +wscloc*gscloc(j,i) & +wliptran*gliptranc(j,i) & +gradafm(j,i) + & +welec*gshieldc(j,i) + & +welec*gshieldc_loc(j,i) + + #else gradc(j,i,icg)=gradbufc(j,i)+welec*gelc(j,i)+ & wel_loc*gel_loc(j,i)+ @@ -723,6 +736,9 @@ c enddo & +wscloc*gscloc(j,i) & +wliptran*gliptranc(j,i) & +gradafm(j,i) + & +welec*gshieldc(j,i) + & +welec*gshieldc_loc(j,i) + #endif gradx(j,i,icg)=wsc*gvdwx(j,i)+wscp*gradx_scp(j,i)+ @@ -731,6 +747,7 @@ c enddo & wsccor*gsccorx(j,i) & +wscloc*gsclocx(j,i) & +wliptran*gliptranx(j,i) + & +welec*gshieldx(j,i) enddo enddo #ifdef DEBUG @@ -3429,7 +3446,9 @@ C do zshift=-1,1 c c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3 c +CTU KURWA do i=iatel_s,iatel_e +C do i=75,75 if (i.le.1) cycle if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1 C changes suggested by Ana to avoid out of bounds @@ -3486,7 +3505,9 @@ c endif c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) num_conti=num_cont_hb(i) +C I TU KURWA do j=ielstart(i),ielend(i) +C do j=16,17 C write (iout,*) i,j if (j.le.1) cycle if (itype(j).eq.ntyp1.or. itype(j+1).eq.ntyp1 @@ -3669,12 +3690,20 @@ c 4/26/02 - AL scaling down 1,4 repulsive VDW interactions el1=fac3*(4.0D0+fac*fac-3.0D0*(cosb*cosb+cosg*cosg)) el2=fac4*fac C MARYSIA - eesij=(el1+el2) +C eesij=(el1+el2) 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 - ees=ees+eesij*fac_shield(i)*fac_shield(j) +C fac_shield(i)=0.4 +C fac_shield(j)=0.6 + el1=el1*fac_shield(i)*fac_shield(j) + el2=el2*fac_shield(i)*fac_shield(j) + eesij=(el1+el2) + ees=ees+eesij else + fac_shield(i)=1.0 + fac_shield(j)=1.0 + eesij=(el1+el2) ees=ees+eesij endif evdw1=evdw1+evdwij*sss @@ -3700,22 +3729,99 @@ C erij(1)=xj*rmij erij(2)=yj*rmij erij(3)=zj*rmij + * * Radial derivatives. First process both termini of the fragment (i,j) * ggg(1)=facel*xj ggg(2)=facel*yj ggg(3)=facel*zj + if ((fac_shield(i).gt.0).and.(fac_shield(j).gt.0).and. + & (shield_mode.gt.0)) then +C print *,i,j + do ilist=1,ishield_list(i) + iresshield=shield_list(ilist,i) + do k=1,3 + rlocshield=grad_shield_side(k,ilist,i)*eesij/fac_shield(i) + gshieldx(k,iresshield)=gshieldx(k,iresshield)+ + & rlocshield + & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i) + gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield +C gshieldc_loc(k,iresshield)=gshieldc_loc(k,iresshield) +C & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i) +C if (iresshield.gt.i) then +C do ishi=i+1,iresshield-1 +C gshieldc(k,ishi)=gshieldc(k,ishi)+rlocshield +C & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i) +C +C enddo +C else +C do ishi=iresshield,i +C gshieldc(k,ishi)=gshieldc(k,ishi)-rlocshield +C & -grad_shield_loc(k,ilist,i)*eesij/fac_shield(i) +C +C enddo +C endif + enddo + enddo + do ilist=1,ishield_list(j) + iresshield=shield_list(ilist,j) + do k=1,3 + rlocshield=grad_shield_side(k,ilist,j)*eesij/fac_shield(j) + gshieldx(k,iresshield)=gshieldx(k,iresshield)+ + & rlocshield + & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j) + gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield + +C & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j) +C gshieldc_loc(k,iresshield)=gshieldc_loc(k,iresshield) +C & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j) +C if (iresshield.gt.j) then +C do ishi=j+1,iresshield-1 +C gshieldc(k,ishi)=gshieldc(k,ishi)+rlocshield +C & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j) +C +C enddo +C else +C do ishi=iresshield,j +C gshieldc(k,ishi)=gshieldc(k,ishi)-rlocshield +C & -grad_shield_loc(k,ilist,j)*eesij/fac_shield(j) +C enddo +C endif + enddo + enddo + + do k=1,3 + gshieldc(k,i)=gshieldc(k,i)+ + & grad_shield(k,i)*eesij/fac_shield(i) + gshieldc(k,j)=gshieldc(k,j)+ + & grad_shield(k,j)*eesij/fac_shield(j) + gshieldc(k,i-1)=gshieldc(k,i-1)+ + & grad_shield(k,i)*eesij/fac_shield(i) + gshieldc(k,j-1)=gshieldc(k,j-1)+ + & grad_shield(k,j)*eesij/fac_shield(j) + + enddo + endif c do k=1,3 c ghalf=0.5D0*ggg(k) c gelc(k,i)=gelc(k,i)+ghalf c gelc(k,j)=gelc(k,j)+ghalf c enddo c 9/28/08 AL Gradient compotents will be summed only at the end +C print *,"before", gelc_long(1,i), gelc_long(1,j) do k=1,3 gelc_long(k,j)=gelc_long(k,j)+ggg(k) +C & +grad_shield(k,j)*eesij/fac_shield(j) gelc_long(k,i)=gelc_long(k,i)-ggg(k) +C & +grad_shield(k,i)*eesij/fac_shield(i) +C gelc_long(k,i-1)=gelc_long(k,i-1) +C & +grad_shield(k,i)*eesij/fac_shield(i) +C gelc_long(k,j-1)=gelc_long(k,j-1) +C & +grad_shield(k,j)*eesij/fac_shield(j) enddo +C print *,"bafter", gelc_long(1,i), gelc_long(1,j) + * * Loop over residues i+1 thru j-1. * @@ -3764,8 +3870,11 @@ C MARYSIA * Radial derivatives. First process both termini of the fragment (i,j) * ggg(1)=fac*xj +C+eesij*grad_shield(1,i)+eesij*grad_shield(1,j) ggg(2)=fac*yj +C+eesij*grad_shield(2,i)+eesij*grad_shield(2,j) ggg(3)=fac*zj +C+eesij*grad_shield(3,i)+eesij*grad_shield(3,j) c do k=1,3 c ghalf=0.5D0*ggg(k) c gelc(k,i)=gelc(k,i)+ghalf @@ -3808,7 +3917,8 @@ c 9/28/08 AL Gradient compotents will be summed only at the end cd print '(2i3,2(3(1pd14.5),3x))',i,j,(dcosb(k),k=1,3), cd & (dcosg(k),k=1,3) do k=1,3 - ggg(k)=ecosb*dcosb(k)+ecosg*dcosg(k) + ggg(k)=(ecosb*dcosb(k)+ecosg*dcosg(k))* + & fac_shield(i)*fac_shield(j) enddo c do k=1,3 c ghalf=0.5D0*ggg(k) @@ -3824,16 +3934,21 @@ cgrad do l=1,3 cgrad gelc(l,k)=gelc(l,k)+ggg(l) cgrad enddo cgrad enddo +C print *,"before22", gelc_long(1,i), gelc_long(1,j) do k=1,3 gelc(k,i)=gelc(k,i) - & +(ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i)) - & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1) + & +((ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i)) + & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)) + & *fac_shield(i)*fac_shield(j) gelc(k,j)=gelc(k,j) - & +(ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j)) - & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1) + & +((ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j)) + & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)) + & *fac_shield(i)*fac_shield(j) gelc_long(k,j)=gelc_long(k,j)+ggg(k) gelc_long(k,i)=gelc_long(k,i)-ggg(k) enddo +C print *,"before33", gelc_long(1,i), gelc_long(1,j) + C MARYSIA c endif !sscale IF (wel_loc.gt.0.0d0 .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 @@ -10530,11 +10645,12 @@ C first for shielding is setting of function of side-chains include 'COMMON.INTERACT' C this is the squar root 77 devided by 81 the epislion in lipid (in protein) double precision div77_81/0.974996043d0/, - &div4_81/0.2222222222d0/ + &div4_81/0.2222222222d0/,sh_frac_dist_grad(3) C the vector between center of side_chain and peptide group double precision pep_side(3),long,side_calf(3), - &pept_group(3) + &pept_group(3),costhet_grad(3),cosphi_grad_long(3), + &cosphi_grad_loc(3),pep_side_norm(3),side_calf_norm(3) C the line belowe needs to be changed for FGPROC>1 do i=1,nres-1 if ((itype(i).eq.ntyp1).and.itype(i+1).eq.ntyp1) cycle @@ -10543,14 +10659,17 @@ Cif there two consequtive dummy atoms there is no peptide group between them C the line below has to be changed for FGPROC>1 VolumeTotal=0.0 do k=1,nres + if ((itype(k).eq.ntyp1).or.(itype(k).eq.10)) cycle dist_pep_side=0.0 dist_side_calf=0.0 do j=1,3 C first lets set vector conecting the ithe side-chain with kth side-chain - pep_side(j)=c(k+nres,j)-(c(i,j)+c(i+1,j))/2.0d0 + pep_side(j)=c(j,k+nres)-(c(j,i)+c(j,i+1))/2.0d0 +C pep_side(j)=2.0d0 C and vector conecting the side-chain with its proper calfa - side_calf(j)=c(k+nres,j)-c(k,j) - pept_group(j)=c(i,j)-c(i+1,j) + side_calf(j)=c(j,k+nres)-c(j,k) +C side_calf(j)=2.0d0 + pept_group(j)=c(j,i)-c(j,i+1) C lets have their lenght dist_pep_side=pep_side(j)**2+dist_pep_side dist_side_calf=dist_side_calf+side_calf(j)**2 @@ -10558,8 +10677,14 @@ C lets have their lenght enddo dist_pep_side=dsqrt(dist_pep_side) dist_pept_group=dsqrt(dist_pept_group) + dist_side_calf=dsqrt(dist_side_calf) + do j=1,3 + pep_side_norm(j)=pep_side(j)/dist_pep_side + side_calf_norm(j)=dist_side_calf + enddo C now sscale fraction sh_frac_dist=-(dist_pep_side-rpp(1,1)-buff_shield)/buff_shield +C print *,buff_shield,"buff" C now sscale if (sh_frac_dist.le.0.0) cycle C If we reach here it means that this side chain reaches the shielding sphere @@ -10567,7 +10692,7 @@ C Lets add him to the list for gradient ishield_list(i)=ishield_list(i)+1 C ishield_list is a list of non 0 side-chain that contribute to factor gradient C this list is essential otherwise problem would be O3 - shield_list(ishield_list(i))=k + shield_list(ishield_list(i),i)=k C Lets have the sscale value if (sh_frac_dist.gt.1.0) then scale_fac_dist=1.0d0 @@ -10577,20 +10702,29 @@ C Lets have the sscale value else scale_fac_dist=-sh_frac_dist*sh_frac_dist & *(2.0*sh_frac_dist-3.0d0) - fac_help_scale=6.0*(scale_fac_dist-scale_fac_dist**2) + fac_help_scale=6.0*(sh_frac_dist-sh_frac_dist**2) & /dist_pep_side/buff_shield*0.5 C remember for the final gradient multiply sh_frac_dist_grad(j) C for side_chain by factor -2 ! do j=1,3 sh_frac_dist_grad(j)=fac_help_scale*pep_side(j) +C print *,"jestem",scale_fac_dist,fac_help_scale, +C & sh_frac_dist_grad(j) enddo endif +C if ((i.eq.3).and.(k.eq.2)) then +C print *,i,sh_frac_dist,dist_pep,fac_help_scale,scale_fac_dist +C & ,"TU" +C endif + C this is what is now we have the distance scaling now volume... short=short_r_sidechain(itype(k)) long=long_r_sidechain(itype(k)) - costhet=1.0d0/dsqrt(1+short**2/dist_pep_side**2) + costhet=1.0d0/dsqrt(1.0+short**2/dist_pep_side**2) C now costhet_grad - costhet_fac=costhet**3*short**2*(-0.5)/dist_pep_side**3 +C costhet=0.0d0 + costhet_fac=costhet**3*short**2*(-0.5)/dist_pep_side**4 +C costhet_fac=0.0d0 do j=1,3 costhet_grad(j)=costhet_fac*pep_side(j) enddo @@ -10602,11 +10736,30 @@ C pep_side0pept_group is vector multiplication do j=1,3 pep_side0pept_group=pep_side0pept_group+pep_side(j)*side_calf(j) enddo - fac_alfa_sin=1.0-(pep_side0pept_group/ - & (dist_pep_side*dist_side_calf))**2 + cosalfa=(pep_side0pept_group/ + & (dist_pep_side*dist_side_calf)) + fac_alfa_sin=1.0-cosalfa**2 fac_alfa_sin=dsqrt(fac_alfa_sin) rkprim=fac_alfa_sin*(long-short)+short - cosphi=1.0d0/dsqrt(1+rkprim**2/dist_pep_side**2) +C now costhet_grad + cosphi=1.0d0/dsqrt(1.0+rkprim**2/dist_pep_side**2) + cosphi_fac=cosphi**3*rkprim**2*(-0.5)/dist_pep_side**4 + + do j=1,3 + cosphi_grad_long(j)=cosphi_fac*pep_side(j) + &+cosphi**3*0.5/dist_pep_side**2*(-rkprim) + &*(long-short)/fac_alfa_sin*cosalfa/ + &((dist_pep_side*dist_side_calf))* + &((side_calf(j))-cosalfa* + &((pep_side(j)/dist_pep_side)*dist_side_calf)) + + cosphi_grad_loc(j)=cosphi**3*0.5/dist_pep_side**2*(-rkprim) + &*(long-short)/fac_alfa_sin*cosalfa + &/((dist_pep_side*dist_side_calf))* + &(pep_side(j)- + &cosalfa*side_calf(j)/dist_side_calf*dist_pep_side) + enddo + VofOverlap=VSolvSphere/2.0d0*(1.0-costhet)*(1.0-cosphi) & /VSolvSphere_div C now the gradient... @@ -10614,27 +10767,29 @@ C grad_shield is gradient of Calfa for peptide groups do j=1,3 grad_shield(j,i)=grad_shield(j,i) C gradient po skalowaniu - & +sh_frac_dist_grad(j)*VofOverlap + & +(sh_frac_dist_grad(j) C gradient po costhet - &+scale_fac_dist*costhet_grad(j) + &-scale_fac_dist*costhet_grad(j)/(1.0-costhet) + &-scale_fac_dist*(cosphi_grad_long(j)) + &/(1.0-cosphi) )*div77_81 + &*VofOverlap C grad_shield_side is Cbeta sidechain gradient grad_shield_side(j,ishield_list(i),i)= - & sh_frac_dist_grad(j)*VofOverlap*2.0d0 - & +scale_fac_dist*costhet_grad(j)*2.0d0 -C grad_shield_side_ca is Calfa sidechain gradient - grad_shield_side_ca(j,ishield_list(i),i)= + & (sh_frac_dist_grad(j)*-2.0d0 + & +scale_fac_dist*costhet_grad(j)*2.0d0/(1.0-costhet) + & +scale_fac_dist*(cosphi_grad_long(j)) + & *2.0d0/(1.0-cosphi)) + & *div77_81*VofOverlap + + grad_shield_loc(j,ishield_list(i),i)= + & scale_fac_dist*cosphi_grad_loc(j) + & *2.0d0/(1.0-cosphi) + & *div77_81*VofOverlap enddo VolumeTotal=VolumeTotal+VofOverlap*scale_fac_dist -C if ((cosphi.le.0.0).or.(costhet.le.0.0)) write(iout,*) "ERROR", -C & cosphi,costhet -C now should be fac_side_grad(k) which will be gradient of factor k which also -C affect the gradient of peptide group i fac_pept_grad(i) and i+1 - write(2,*) "myvolume",VofOverlap,VSolvSphere_div,VolumeTotal - enddo -C write(2,*) "TOTAL VOLUME",i,VolumeTotal -C the scaling factor of the shielding effect + enddo fac_shield(i)=VolumeTotal*div77_81+div4_81 - write(2,*) "TOTAL VOLUME",i,VolumeTotal,fac_shield(i) +C write(2,*) "TOTAL VOLUME",i,VolumeTotal,fac_shield(i) enddo return end diff --git a/source/unres/src_MD-M/gradient_p.F b/source/unres/src_MD-M/gradient_p.F index 9961837..edbf1bc 100644 --- a/source/unres/src_MD-M/gradient_p.F +++ b/source/unres/src_MD-M/gradient_p.F @@ -346,6 +346,7 @@ C------------------------------------------------------------------------- include 'COMMON.VAR' include 'COMMON.MD' include 'COMMON.SCCOR' + include 'COMMON.SHIELD' C C Initialize Cartesian-coordinate gradient C @@ -357,6 +358,9 @@ C gvdwc_scp(j,i)=0.0D0 gvdwc_scpp(j,i)=0.0d0 gelc (j,i)=0.0D0 + gshieldx(j,i)=0.0d0 + gshieldc(j,i)=0.0d0 + gshieldc_loc(j,i)=0.0d0 gelc_long(j,i)=0.0D0 gradb(j,i)=0.0d0 gradbx(j,i)=0.0d0 @@ -388,8 +392,12 @@ C C grad_shield_side is Cbeta sidechain gradient do kk=1,maxshieldlist grad_shield_side(j,kk,i)=0.0d0 + grad_shield_loc(j,kk,i)=0.0d0 + C grad_shield_side_ca is Calfa sidechain gradient - grad_shield_side_ca(j,kk,i)=0.0d0 + + +C grad_shield_side_ca(j,kk,i)=0.0d0 enddo do intertyp=1,3 gloc_sc(intertyp,i,icg)=0.0d0 diff --git a/source/unres/src_MD-M/parmread.F b/source/unres/src_MD-M/parmread.F index 574552c..c881425 100644 --- a/source/unres/src_MD-M/parmread.F +++ b/source/unres/src_MD-M/parmread.F @@ -949,14 +949,14 @@ c B2(1,-i) =b(2) c B2(2,-i) =-b(4) B1tilde(1,i) = b(3,i) B1tilde(2,i) =-b(5,i) - B1tilde(1,-i) =-b(3,i) - B1tilde(2,-i) =b(5,i) +C B1tilde(1,-i) =-b(3,i) +C B1tilde(2,-i) =b(5,i) b1tilde(1,i)=0.0d0 b1tilde(2,i)=0.0d0 B2(1,i) = b(2,i) B2(2,i) = b(4,i) - B2(1,-i) =b(2,i) - B2(2,-i) =-b(4,i) +C B2(1,-i) =b(2,i) +C B2(2,-i) =-b(4,i) c b2(1,i)=0.0d0 c b2(2,i)=0.0d0 diff --git a/source/unres/src_MD-M/readrtns_CSA.F b/source/unres/src_MD-M/readrtns_CSA.F index 4d70e89..bdc7c9d 100644 --- a/source/unres/src_MD-M/readrtns_CSA.F +++ b/source/unres/src_MD-M/readrtns_CSA.F @@ -274,6 +274,7 @@ C long axis of side chain long_r_sidechain(i)=vbldsc0(1,i) short_r_sidechain(i)=sigma0(i) enddo + buff_shield=1.0d0 endif if (me.eq.king .or. .not.out1file ) & write (iout,*) "DISTCHAINMAX",distchainmax