X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fenergy.F90;h=63a76fb11c330f5f472172d0987d47bf6cdc4482;hb=f7883bb11fc7df1ed921187999982670d26e6fa5;hp=376d3181932b1983f6001b5b45dcda773fb58e44;hpb=94d31d7bcffc5412d9a88a0aef46b68349b60cf3;p=unres4.git diff --git a/source/unres/energy.F90 b/source/unres/energy.F90 index 376d318..63a76fb 100644 --- a/source/unres/energy.F90 +++ b/source/unres/energy.F90 @@ -2007,7 +2007,7 @@ zj=c(3,nres+j) call to_box(xj,yj,zj) call lipid_layer(xj,yj,zj,sslipj,ssgradlipj) - write (iout,*) "KWA2", itypi,itypj +! write (iout,*) "KWA2", itypi,itypj aa=aa_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 bb=bb_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & @@ -5555,7 +5555,7 @@ iabs(itype(jjj,1)).eq.1) then call ssbond_ene(iii,jjj,eij) ehpb=ehpb+2*eij -!d write (iout,*) "eij",eij +! write (iout,*) "eij",eij,iii,jjj endif else if (ii.gt.nres .and. jj.gt.nres) then !c Restraints from contact prediction @@ -5694,10 +5694,13 @@ itypj=iabs(itype(j,1)) ! dscj_inv=dsc_inv(itypj) dscj_inv=vbld_inv(nres+j) - xj=c(1,nres+j)-xi - yj=c(2,nres+j)-yi - zj=c(3,nres+j)-zi + xj=c(1,nres+j) + yj=c(2,nres+j) + zj=c(3,nres+j) call to_box(xj,yj,zj) + xj=boxshift(xj-xi,boxxsize) + yj=boxshift(yj-yi,boxysize) + zj=boxshift(zj-zi,boxzsize) dxj=dc_norm(1,nres+j) dyj=dc_norm(2,nres+j) dzj=dc_norm(3,nres+j) @@ -5722,9 +5725,9 @@ eij=akcm*deltad*deltad+akth*(deltat1*deltat1+deltat2*deltat2) & +akct*deltad*deltat12 & +v1ss*cosphi+v2ss*cosphi*cosphi+v3ss*cosphi*cosphi*cosphi+ebr -! write(iout,*) i,j,"rij",rij,"d0cm",d0cm," akcm",akcm," akth",akth, -! & " akct",akct," deltad",deltad," deltat",deltat1,deltat2, -! & " deltat12",deltat12," eij",eij +! write(iout,*) i,j,"rij",rij,"d0cm",d0cm," akcm",akcm," akth",akth, & +! " akct",akct," deltad",deltad," deltat",deltat1,deltat2, & +! " deltat12",deltat12," eij",eij ed=2*akcm*deltad+akct*deltat12 pom1=akct*deltad pom2=v1ss+2*v2ss*cosphi+3*v3ss*cosphi*cosphi @@ -13412,6 +13415,10 @@ yj=c(2,nres+j)-yi zj=c(3,nres+j)-zi call to_box(xj,yj,zj) + xj=boxshift(xj-xi,boxxsize) + yj=boxshift(yj-yi,boxysize) + zj=boxshift(zj-zi,boxzsize) + rrij=1.0D0/(xj*xj+yj*yj+zj*zj) fac_augm=rrij**expon e_augm=augm(itypi,itypj)*fac_augm @@ -22271,6 +22278,10 @@ chip1=chip(itypi) real(kind=8) :: facd4, adler, Fgb, facd3 integer troll,jj,istate real (kind=8) :: dcosom1(3),dcosom2(3) + real(kind=8) ::locbox(3) + locbox(1)=boxxsize + locbox(2)=boxysize + locbox(3)=boxzsize evdw=0.0D0 if (nres_molec(5).eq.0) return @@ -22313,6 +22324,8 @@ chip1=chip(itypi) zj=c(3,j) call to_box(xj,yj,zj) +! write(iout,*) "xi,yi,zi,xj,yj,zj", xi,yi,zi,xj,yj,zj + ! call lipid_layer(xj,yj,zj,sslipj,ssgradlipj) ! aa=aa_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & ! +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 @@ -22321,6 +22334,7 @@ chip1=chip(itypi) xj=boxshift(xj-xi,boxxsize) yj=boxshift(yj-yi,boxysize) zj=boxshift(zj-zi,boxzsize) +! write(iout,*) "xj,yj,zj", xj,yj,zj,boxxsize ! dxj = dc_norm( 1, nres+j ) ! dyj = dc_norm( 2, nres+j ) @@ -22364,11 +22378,13 @@ chip1=chip(itypi) ctail(k,1)=c(k,i+nres) ctail(k,2)=c(k,j) END DO + call to_box(ctail(1,1),ctail(2,1),ctail(3,1)) + call to_box(ctail(1,2),ctail(2,2),ctail(3,2)) !c! tail distances will be themselves usefull elswhere !c1 (in Gcav, for example) - Rtail_distance(1) = ctail( 1, 2 ) - ctail( 1,1 ) - Rtail_distance(2) = ctail( 2, 2 ) - ctail( 2,1 ) - Rtail_distance(3) = ctail( 3, 2 ) - ctail( 3,1 ) + do k=1,3 + Rtail_distance(k) = boxshift(ctail(k,2) - ctail(k,1),locbox(k)) + enddo Rtail = dsqrt( & (Rtail_distance(1)*Rtail_distance(1)) & + (Rtail_distance(2)*Rtail_distance(2)) & @@ -22383,10 +22399,15 @@ chip1=chip(itypi) ! see unres publications for very informative images chead(k,1) = c(k, i+nres) + d1 * dc_norm(k, i+nres) chead(k,2) = c(k, j) + enddo + call to_box(chead(1,1),chead(2,1),chead(3,1)) + call to_box(chead(1,2),chead(2,2),chead(3,2)) + ! distance ! Rsc_distance(k) = dabs(c(k, i+nres) - c(k, j+nres)) -! Rsc(k) = Rsc_distance(k) * Rsc_distance(k) - Rhead_distance(k) = chead(k,2) - chead(k,1) +! Rsc(k) = Rsc_distance(k) * Rsc_distance(k) + do k=1,3 + Rhead_distance(k) = boxshift(chead(k,2) - chead(k,1),locbox(k)) END DO ! pitagoras (root of sum of squares) Rhead = dsqrt( & @@ -22408,6 +22429,7 @@ chip1=chip(itypi) dPOLdOM1 = 0.0d0 dPOLdOM2 = 0.0d0 Fcav = 0.0d0 + Fisocav=0.0d0 dFdR = 0.0d0 dCAVdOM1 = 0.0d0 dCAVdOM2 = 0.0d0 @@ -22635,6 +22657,10 @@ chip1=chip(itypi) yj=c(2,j) zj=c(3,j) call to_box(xj,yj,zj) + xj=boxshift(xj-xi,boxxsize) + yj=boxshift(yj-yi,boxysize) + zj=boxshift(zj-zi,boxzsize) + dist_init=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2 dxj = 0.0d0! dc_norm( 1, nres+j ) @@ -22679,11 +22705,16 @@ chip1=chip(itypi) ctail(k,1)=(c(k,i)+c(k,i+1))/2.0 ctail(k,2)=c(k,j) END DO + call to_box(ctail(1,1),ctail(2,1),ctail(3,1)) + call to_box(ctail(1,2),ctail(2,2),ctail(3,2)) +!c! tail distances will be themselves usefull elswhere +!c1 (in Gcav, for example) + do k=1,3 + Rtail_distance(k) = boxshift(ctail(k,2) - ctail(k,1),locbox(k)) + enddo + !c! tail distances will be themselves usefull elswhere !c1 (in Gcav, for example) - Rtail_distance(1) = ctail( 1, 2 ) - ctail( 1,1 ) - Rtail_distance(2) = ctail( 2, 2 ) - ctail( 2,1 ) - Rtail_distance(3) = ctail( 3, 2 ) - ctail( 3,1 ) Rtail = dsqrt( & (Rtail_distance(1)*Rtail_distance(1)) & + (Rtail_distance(2)*Rtail_distance(2)) & @@ -22700,11 +22731,20 @@ chip1=chip(itypi) ! see unres publications for very informative images chead(k,1) = (c(k, i)+c(k,i+1))/2.0 + d1 * dc_norm(k, i) chead(k,2) = c(k, j) + ENDDO ! distance ! Rsc_distance(k) = dabs(c(k, i+nres) - c(k, j+nres)) ! Rsc(k) = Rsc_distance(k) * Rsc_distance(k) - Rhead_distance(k) = chead(k,2) - chead(k,1) + call to_box(chead(1,1),chead(2,1),chead(3,1)) + call to_box(chead(1,2),chead(2,2),chead(3,2)) + +! distance +! Rsc_distance(k) = dabs(c(k, i+nres) - c(k, j+nres)) +! Rsc(k) = Rsc_distance(k) * Rsc_distance(k) + do k=1,3 + Rhead_distance(k) = boxshift(chead(k,2) - chead(k,1),locbox(k)) END DO + ! pitagoras (root of sum of squares) Rhead = dsqrt( & (Rhead_distance(1)*Rhead_distance(1)) & @@ -23453,6 +23493,7 @@ chip1=chip(itypi) yj=c(2,j) zj=c(3,j) call to_box(xj,yj,zj) +! write(iout,*) "xi,yi,zi,xj,yj,zj", xi,yi,zi,xj,yj,zj ! call lipid_layer(xj,yj,zj,sslipj,ssgradlipj) ! aa=aa_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & ! +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 @@ -23461,7 +23502,7 @@ chip1=chip(itypi) xj=boxshift(xj-xi,boxxsize) yj=boxshift(yj-yi,boxysize) zj=boxshift(zj-zi,boxzsize) - +! write(iout,*) 'after shift', xj,yj,zj dist_init=xj**2+yj**2+zj**2 itypi=itype(i,2) @@ -25240,28 +25281,7 @@ chip1=chip(itypi) zi=c(3,nres+i) call to_box(xi,yi,zi) call lipid_layer(xi,yi,zi,sslipi,ssgradlipi) - if ((zi.gt.bordlipbot) & - .and.(zi.lt.bordliptop)) then -!C the energy transfer exist - if (zi.lt.buflipbot) then -!C what fraction I am in - fracinbuf=1.0d0- & - ((zi-bordlipbot)/lipbufthick) -!C lipbufthick is thickenes of lipid buffore - sslipi=sscalelip(fracinbuf) - ssgradlipi=-sscagradlip(fracinbuf)/lipbufthick - elseif (zi.gt.bufliptop) then - fracinbuf=1.0d0-((bordliptop-zi)/lipbufthick) - sslipi=sscalelip(fracinbuf) - ssgradlipi=sscagradlip(fracinbuf)/lipbufthick - else - sslipi=1.0d0 - ssgradlipi=0.0 - endif - else - sslipi=0.0d0 - ssgradlipi=0.0 - endif +! endif ! print *, sslipi,ssgradlipi dxi=dc_norm(1,nres+i) dyi=dc_norm(2,nres+i) @@ -25318,7 +25338,7 @@ chip1=chip(itypi) zj=c(3,j+nres) call to_box(xj,yj,zj) call lipid_layer(xj,yj,zj,sslipj,ssgradlipj) - write(iout,*) "KRUWA", i,j +! write(iout,*) "KRUWA", i,j aa=aa_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 & +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0 bb=bb_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0 &