From cc980c20ed89622e129591e212a044280de43b9f Mon Sep 17 00:00:00 2001 From: Adam Kazimierz Sieradzan Date: Tue, 29 May 2012 05:54:26 -0400 Subject: [PATCH] sadas --- source/unres/src_MD/COMMON.SCCOR | 13 +- source/unres/src_MD/COMMON.VAR | 1 + source/unres/src_MD/checkder_p.F | 15 +- source/unres/src_MD/energy_p_new_barrier.F | 36 +++-- source/unres/src_MD/int_to_cart.f | 102 ++++++++++++ source/unres/src_MD/intcartderiv.F | 241 ++++++++++++++++++++++++++++ source/unres/src_MD/parmread.F | 52 ++++-- 7 files changed, 432 insertions(+), 28 deletions(-) diff --git a/source/unres/src_MD/COMMON.SCCOR b/source/unres/src_MD/COMMON.SCCOR index 5217de7..1c725d7 100644 --- a/source/unres/src_MD/COMMON.SCCOR +++ b/source/unres/src_MD/COMMON.SCCOR @@ -1,6 +1,11 @@ -C Parameters of the SCCOR term - double precision v1sccor,v2sccor - integer nterm_sccor + Parameters of the SCCOR term + double precision v1sccor,v2sccor,vlor1sccor, + & vlor2sccor,vlor3sccor + integer nterm_sccor,isccortyp,nsccortyp,nlor_sccor common/torsion/v1sccor(maxterm_sccor,20,20), & v2sccor(maxterm_sccor,20,20), - & nterm_sccor + & nterm_sccor(ntyp,ntyp),isccortyp(ntyp),nsccortyp, + & nlor_sccor(ntyp,ntyp),vlor1sccor(maxterm_sccor,20,20), + & vlor2sccor(maxterm_sccor,20,20), + & vlor3sccor(maxterm_sccor,20,20) + diff --git a/source/unres/src_MD/COMMON.VAR b/source/unres/src_MD/COMMON.VAR index 71158b8..c2a0a25 100644 --- a/source/unres/src_MD/COMMON.VAR +++ b/source/unres/src_MD/COMMON.VAR @@ -5,6 +5,7 @@ C Store the geometric variables in the following COMMON block. & thetaref,phiref,costtab,sinttab,cost2tab,sint2tab, & xxtab,yytab,zztab,xxref,yyref,zzref common /var/ theta(maxres),phi(maxres),alph(maxres),omeg(maxres), + & omicron(2,maxres),tauangle(3,maxres) & vbld(2*maxres),thetaref(maxres),phiref(maxres), & costtab(maxres), sinttab(maxres), cost2tab(maxres), & sint2tab(maxres),xxtab(maxres),yytab(maxres), diff --git a/source/unres/src_MD/checkder_p.F b/source/unres/src_MD/checkder_p.F index 26854e6..ec66acd 100644 --- a/source/unres/src_MD/checkder_p.F +++ b/source/unres/src_MD/checkder_p.F @@ -513,7 +513,20 @@ c------------------------------------------------------------------------- & +(c(j,i+1)-c(j,i))/dnorm2) enddo be=0.0D0 - if (i.gt.2) phi(i+1)=beta(i-2,i-1,i,i+1) + if (i.gt.2) then + phi(i+1)=beta(i-2,i-1,i,i+1) + if (itype(i).ne.10).and.(itype(i-1).ne.10) then + tauangle(3,i+1)=beta(i+nres-1,i-1,i,i+nres) + endif + if (itype(i-1).ne.10) then + tauangle(1,i+1)=beta(i-1+nres,i-1,i,i+1) + omicron(1,i)=alpha(i-2,i-1,i-1+nres) + omicron(2,i)=alpha(i-1+nres,i-1,i) + endif + if (itype(i).ne.10) then + tauangle(2,i+1)=beta(i-2,i-1,i,i+nres) + endif + endif omeg(i)=beta(nres+i,i,maxres2,i+1) alph(i)=alpha(nres+i,i,maxres2) theta(i+1)=alpha(i-1,i,i+1) diff --git a/source/unres/src_MD/energy_p_new_barrier.F b/source/unres/src_MD/energy_p_new_barrier.F index e3f2c36..23449da 100644 --- a/source/unres/src_MD/energy_p_new_barrier.F +++ b/source/unres/src_MD/energy_p_new_barrier.F @@ -5801,7 +5801,6 @@ c lprn=.true. phii1=phi(i+1) gloci1=0.0D0 gloci2=0.0D0 -C Regular cosine and sine terms do j=1,ntermd_1(itori,itori1,itori2) v1cij=v1c(1,j,itori,itori1,itori2) v1sij=v1s(1,j,itori,itori1,itori2) @@ -5870,24 +5869,43 @@ c write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor esccor=0.0D0 do i=iphi_start,iphi_end esccor_ii=0.0D0 - itori=itype(i-2) - itori1=itype(i-1) + isccori=isccortyp(itype(i-2)) + isccori1=isccortyp(itype(i-1)) phii=phi(i) +cccc Added 9 May 2012 +cc Tauangle is torsional engle depending on the value of first digit +c(see comment below) +cc Omicron is flat angle depending on the value of first digit +c(see comment below) + gloci=0.0D0 - do j=1,nterm_sccor - v1ij=v1sccor(j,itori,itori1) - v2ij=v2sccor(j,itori,itori1) - cosphi=dcos(j*phii) - sinphi=dsin(j*phii) + do intertyp=1,3 +cc Added 09 May 2012 (Adasko) +cc Intertyp means interaction type of backbone mainchain correlation: +c 1 = SC...Ca...Ca...Ca +c 2 = Ca...Ca...Ca...SC +c 3 = SC...Ca...Ca...SC + if (((intertyp.eq.3).and.(itype(i-2).eq.10).or. + & (itype(i-1).eq.10)) + & .or. ((intertyp.eq.1).and.(itype(i-2).ne.10)) + & .or. ((intertyp.eq.2).and.(itype(i-1).ne.10))) cycle + do j=1,nterm_sccor(isccori,isccori1) + v1ij=v1sccor(j,intertyp,isccori,isccori1) + v2ij=v2sccor(j,intertyp,isccori,isccori1) + cosphi=dcos(j*tauangle(intertyp,i)) + sinphi=dsin(j*tauangle(intertyp,i)) esccor=esccor+v1ij*cosphi+v2ij*sinphi gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi) enddo + gloc_sc(intertyp,i-3,icg)=gloc_sc(i-3,icg)+wtor*gloci if (lprn) & write (iout,'(2(a3,2x,i3,2x),2i3,6f8.3/26x,6f8.3/)') & restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,itori,itori1, - & (v1sccor(j,itori,itori1),j=1,6),(v2sccor(j,itori,itori1),j=1,6) + & (v1sccor(j,intertyp,itori,itori1),j=1,6) + & ,(v2sccor(j,intertyp,itori,itori1),j=1,6) gsccor_loc(i-3)=gsccor_loc(i-3)+gloci enddo + enddo return end c---------------------------------------------------------------------------- diff --git a/source/unres/src_MD/int_to_cart.f b/source/unres/src_MD/int_to_cart.f index 0528af7..2564372 100644 --- a/source/unres/src_MD/int_to_cart.f +++ b/source/unres/src_MD/int_to_cart.f @@ -113,6 +113,108 @@ c The side-chain vector derivatives enddo endif enddo +C INTERTYP=1 SC...Ca...Ca...Ca +c calculating dE/ddc1 + intertyp=1 + if (nres.lt.2) return + if ((nres.lt.3).and.(itype(1).eq.1)) return + if (itype(1).ne.10) then + do j=1,3 + gxcart(j,1)=gxcart(j,1)+gloc_sc(intertyp,1,icg)* + & dtauangle(j,1,1,4) +c As potetnial DO NOT +c & +gloc_sc(intertyp,nres-2,icg)*dtheta(j,1,3) + enddo + endif +c Calculating the remainder of dE/ddc2 + do j=1,3 + gcart(j,2)=gcart(j,2)+gloc(1,icg)*dphi(j,2,4)+ + & gloc(nres-2,icg)*dtheta(j,2,3)+gloc(nres-1,icg)*dtheta(j,1,4) + if(itype(2).ne.10) then + gcart(j,2)=gcart(j,2)+gloc(ialph(2,1),icg)*dalpha(j,2,2)+ + & gloc(ialph(2,1)+nside,icg)*domega(j,2,2) + endif + if(itype(3).ne.10) then + gcart(j,2)=gcart(j,2)+gloc(ialph(3,1),icg)*dalpha(j,1,3)+ + & gloc(ialph(3,1)+nside,icg)*domega(j,1,3) + endif + if(nres.gt.4) then + gcart(j,2)=gcart(j,2)+gloc(2,icg)*dphi(j,1,5) + endif + enddo +c If there are only five residues + if(nres.eq.5) then + do j=1,3 + gcart(j,3)=gcart(j,3)+gloc(1,icg)*dphi(j,3,4)+gloc(2,icg)* + & dphi(j,2,5)+gloc(nres-1,icg)*dtheta(j,2,4)+gloc(nres,icg)* + & dtheta(j,1,5) + if(itype(3).ne.10) then + gcart(j,3)=gcart(j,3)+gloc(ialph(3,1),icg)* + & dalpha(j,2,3)+gloc(ialph(3,1)+nside,icg)*domega(j,2,3) + endif + if(itype(4).ne.10) then + gcart(j,3)=gcart(j,3)+gloc(ialph(4,1),icg)* + & dalpha(j,1,4)+gloc(ialph(4,1)+nside,icg)*domega(j,1,4) + endif + enddo + endif +c If there are more than five residues + if(nres.gt.5) then + do i=3,nres-3 + do j=1,3 + gcart(j,i)=gcart(j,i)+gloc(i-2,icg)*dphi(j,3,i+1) + & +gloc(i-1,icg)*dphi(j,2,i+2)+ + & gloc(i,icg)*dphi(j,1,i+3)+gloc(nres+i-4,icg)*dtheta(j,2,i+1)+ + & gloc(nres+i-3,icg)*dtheta(j,1,i+2) + if(itype(i).ne.10) then + gcart(j,i)=gcart(j,i)+gloc(ialph(i,1),icg)*dalpha(j,2,i)+ + & gloc(ialph(i,1)+nside,icg)*domega(j,2,i) + endif + if(itype(i+1).ne.10) then + gcart(j,i)=gcart(j,i)+gloc(ialph(i+1,1),icg)*dalpha(j,1,i+1) + & +gloc(ialph(i+1,1)+nside,icg)*domega(j,1,i+1) + endif + enddo + enddo + endif +c Setting dE/ddnres-2 + if(nres.gt.5) then + do j=1,3 + gcart(j,nres-2)=gcart(j,nres-2)+gloc(nres-4,icg)* + & dphi(j,3,nres-1)+gloc(nres-3,icg)*dphi(j,2,nres) + & +gloc(2*nres-6,icg)* + & dtheta(j,2,nres-1)+gloc(2*nres-5,icg)*dtheta(j,1,nres) + if(itype(nres-2).ne.10) then + gcart(j,nres-2)=gcart(j,nres-2)+gloc(ialph(nres-2,1),icg)* + & dalpha(j,2,nres-2)+gloc(ialph(nres-2,1)+nside,icg)* + & domega(j,2,nres-2) + endif + if(itype(nres-1).ne.10) then + gcart(j,nres-2)=gcart(j,nres-2)+gloc(ialph(nres-1,1),icg)* + & dalpha(j,1,nres-1)+gloc(ialph(nres-1,1)+nside,icg)* + & domega(j,1,nres-1) + endif + enddo + endif +c Settind dE/ddnres-1 + do j=1,3 + gcart(j,nres-1)=gcart(j,nres-1)+gloc(nres-3,icg)*dphi(j,3,nres)+ + & gloc(2*nres-5,icg)*dtheta(j,2,nres) + if(itype(nres-1).ne.10) then + gcart(j,nres-1)=gcart(j,nres-1)+gloc(ialph(nres-1,1),icg)* + & dalpha(j,2,nres-1)+gloc(ialph(nres-1,1)+nside,icg)* + & domega(j,2,nres-1) + endif + enddo +c The side-chain vector derivatives + do i=2,nres-1 + if(itype(i).ne.10) then + do j=1,3 + gxcart(j,i)=gxcart(j,i)+gloc(ialph(i,1),icg)*dalpha(j,3,i) + & +gloc(ialph(i,1)+nside,icg)*domega(j,3,i) + enddo + endif + enddo return end diff --git a/source/unres/src_MD/intcartderiv.F b/source/unres/src_MD/intcartderiv.F index 5fea875..c7cf95c 100644 --- a/source/unres/src_MD/intcartderiv.F +++ b/source/unres/src_MD/intcartderiv.F @@ -44,6 +44,43 @@ c We need dtheta(:,:,i-1) to compute dphi(:,:,i) dtheta(j,2,i)=-1/sint*dcostheta(j,2,i) enddo enddo + +#if defined(MPI) && defined(PARINTDER) +c We need dtheta(:,:,i-1) to compute dphi(:,:,i) + do i=max0(ithet_start-1,3),ithet_end +#else + do i=3,nres +#endif + if (itype(i-1).ne.10) then + cost1=dcos(omicron(1,i)) + sint1=sqrt(1-cost1*cost1) + cost2=dcos(omicron(2,i)) + sint2=sqrt(1-cost2*cost2) + do j=1,3 +CC Calculate derivative over first omicron (Cai-2,Cai-1,SCi-1) + dcosomicron(j,1,1,i)=-(-dc_norm(j,i-1+nres)+ + & cost1*dc_norm(j,i-2))/ + & vbld(i-1) + domicron(j,1,1,i)=-1/sint1*dcosomicron(j,1,1,i) + dcosomicron(j,1,2,i)=-(dc_norm(j,i-2) + & +cost1*(-dc_norm(j,i-1+nres)))/ + & vbld(i+nres) + domicron(j,1,2,i)=-1/sint1*dcosomicron(j,1,2,i) +CC Calculate derivative over second omicron Sci-1,Cai-1 Cai +CC Looks messy but better than if in loop + dcosomicron(j,2,1,i)=-(-dc_norm(j,i-1+nres) + & +cost2*dc_norm(j,i-1))/ + & vbld(i) + domicron(j,2,1,i)=-1/sint2*dcosomicron(j,2,1,i) + dcosomicron(j,2,2,i)=-(dc_norm(j,i-1) + & +cost1*(-dc_norm(j,i-1+nres)))/ + & vbld(i+nres) + domicron(j,2,2,i)=-1/sint2*dcosomicron(j,2,2,i) + enddo + endif + enddo + + c Derivatives of phi: c If phi is 0 or 180 degrees, then the formulas @@ -109,6 +146,210 @@ c Obtaining the gamma derivatives from cosine derivative enddo endif enddo + +Calculate derivative of Tauangle +#ifdef PARINTDER + do i=iphi1_start,iphi1_end +#else + do i=4,nres +#endif +cc INTERTYP=1 SC...Ca...Ca..Ca +c the conventional case + sint=dsin(theta(i)) + sint1=dsin(omicron(2,i-1)) + sing=dsin(tauangle(1,i)) + cost=dcos(theta(i)) + cost1=dcos(omicron(2,i-1)) + cosg=dcos(tauangle(1,i)) + do j=1,3 + dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres) + enddo + scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1)) + fac0=1.0d0/(sint1*sint) + fac1=cost*fac0 + fac2=cost1*fac0 + fac3=cosg*cost1/(sint1*sint1) + fac4=cosg*cost/(sint*sint) +c Obtaining the gamma derivatives from sine derivative + if (tauangle(1,i).gt.-pi4.and.tauangle(1,i).le.pi4.or. + & tauangle(1,i).gt.pi34.and.tauangle(1,i).le.pi.or. + & tauangle(1,i).gt.-pi.and.tauangle(1,i).le.-pi34) then + do j=1,3 + dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres) + enddo + call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1) + call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-1),vp2) + call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3) + do j=1,3 + ctgt=cost/sint + ctgt1=cost1/sint1 + cosg_inv=1.0d0/cosg + dsintau(j,1,1,i)=-sing*ctgt1*domicron(j,2,1,i-1) + &-(fac0*vp1(j)+sing*(-dc_norm(j,i-2+nres))) + & *vbld_inv(i-2+nres) + dtauangle(j,1,1,i)=cosg_inv*dsintau(j,1,1,i) + dsintau(j,1,2,i)= + & -sing*(ctgt1*domicron(j,2,2,i-1)+ctgt*dtheta(j,1,i)) + & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1) + dtauangle(j,1,2,i)=cosg_inv*dsintau(j,1,2,i) +c Bug fixed 3/24/05 (AL) + dsintau(j,1,3,i)=-sing*ctgt*dtheta(j,2,i) + & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i) +c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1) + dtauangle(j,1,3,i)=cosg_inv*dsinphi(j,1,3,i) + enddo +c Obtaining the gamma derivatives from cosine derivative + else + do j=1,3 + dcostau(j,1,1,i)=fac1*dcosomicron(j,2,1,i-1)+fac3* + & dcosomicron(j,2,1,i-1)-fac0*(dc_norm(j,i-1)-scalp* + & (-dc_norm(j,i-2+nres)))/vbld(i-2+nres) + dtauangle(j,1,1,i)=-1/sing*dcostau(j,1,1,i) + dcostau(j,1,2,i)=fac1*dcosomicron(j,2,2,i-1)+fac2* + & dcostheta(j,1,i)+fac3*dcosomicron(j,2,2,i-1)+fac4* + & dcostheta(j,1,i) + dtauangle(j,1,2,i)=-1/sing*dcostau(j,1,2,i) + dcostau(j,1,3,i)=fac2*dcostheta(j,2,i)+fac4* + & dcostheta(j,2,i)-fac0*(-dc_norm(j,i-2+nres)-scalp* + & dc_norm(j,i-1))/vbld(i) + dtauangle(j,1,3,i)=-1/sing*dcostau(j,1,3,i) + enddo + endif + enddo +CC Second case Ca...Ca...Ca...SC +#ifdef PARINTDER + do i=iphi1_start,iphi1_end +#else + do i=4,nres +#endif +c the conventional case + sint=dsin(omicron(1,i)) + sint1=dsin(theta(i-1)) + sing=dsin(tauangle(2,i) + cost=dcos(omicron(1,i) + cost1=dcos(theta(i-1)) + cosg=dcos(tauangle(2,i)) + do j=1,3 + dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres) + enddo + scalp=scalar(dc_norm(1,i-3),dc_norm2(1,i-1+nres)) + fac0=1.0d0/(sint1*sint) + fac1=cost*fac0 + fac2=cost1*fac0 + fac3=cosg*cost1/(sint1*sint1) + fac4=cosg*cost/(sint*sint) +c Obtaining the gamma derivatives from sine derivative + if (tauangle(2,i).gt.-pi4.and.tauangle(2,i).le.pi4.or. + & tauangle(2,i).gt.pi34.and.tauangle(2,i).le.pi.or. + & tauangle(2,i).gt.-pi.and.tauangle(2,i).le.-pi34) then + call vecpr(dc_norm2(1,i-1+nres),dc_norm(1,i-2),vp1) + call vecpr(dc_norm(1,i-3),dc_norm2(1,i-1+nres),vp2) + call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3) + do j=1,3 + ctgt=cost/sint + ctgt1=cost1/sint1 + cosg_inv=1.0d0/cosg + dsintau(j,2,1,i)=-sing*ctgt1*dtheta(j,1,i-1) + & -(fac0*vp1(j)+sing*dc_norm(j,i-3))*vbld_inv(i-2) + dtauangle(j,2,1,i)=cosg_inv*dsintau(j,2,1,i) + dsintau(j,2,2,i)= + & -sing*(ctgt1*dtheta(j,2,i-1)+ctgt*domicron(j,1,1,i)) + & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1) + dphi(j,2,i)=cosg_inv*dsintau(j,2,2,i) +c Bug fixed 3/24/05 (AL) + dsintau(j,2,3,i)=-sing*ctgt*domicron(j,1,2,i) + & +(fac0*vp3(j)-sing*dc_norm2(j,i-1+nres))*vbld_inv(i) +c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1) + dtauangle(j,2,3,i)=cosg_inv*dsinphi(j,2,3,i) + enddo +c Obtaining the gamma derivatives from cosine derivative + else + do j=1,3 + dcostau(j,2,1,i)=fac1*dcostheta(j,1,i-1)+fac3* + & dcostheta(j,1,i-1)-fac0*(dc_norm2(j,i-1+nres)-scalp* + & dc_norm(j,i-3))/vbld(i-2) + dtauangle(j,2,1,i)=-1/sing*dcostau(j,2,1,i) + dcostau(j,2,2,i)=fac1*dcostheta(j,2,i-1)+fac2* + & dcosomicron(j,1,1,i)+fac3*dcostheta(j,2,i-1)+fac4* + & dcosomicron(j,1,1,i) + dtauanlge(j,2,2,i)=-1/sing*dcostau(j,2,2,i) + dcostau(j,2,3,i)=fac2*dcosomicron(j,1,2,i)+fac4* + & dcostheta(j,1,2,i)-fac0*(dc_norm(j,i-3)-scalp* + & dc_norm2(j,i-1+nres))/vbld(i-1+nres) + dtauanlge(j,2,3,i)=-1/sing*dcosphi(j,3,i) + enddo + endif + enddo + + +CCC third case SC...Ca...Ca...SC +#ifdef PARINTDER + + do i=iphi1_start,iphi1_end +#else + do i=4,nres +#endif +c the conventional case + sint=dsin(omicron(1,i)) + sint1=dsin(omicron(2,i-1)) + sing=dsin(tauangle(3,i)) + cost=dcos(omicron(1,i)) + cost1=dcos(omicron(2,i-1)) + cosg=dcos(tauangle(3,i)) + do j=1,3 + dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres) + dc_norm2(j,i-1+nres)=-dc_norm(j,i-2+nres) + enddo + scalp=scalar(dc_norm2(1,i-2+nres),dc_norm2(1,i-1+nres)) + fac0=1.0d0/(sint1*sint) + fac1=cost*fac0 + fac2=cost1*fac0 + fac3=cosg*cost1/(sint1*sint1) + fac4=cosg*cost/(sint*sint) +c Obtaining the gamma derivatives from sine derivative + if (tauangle(3,i).gt.-pi4.and.tauangle(3,i).le.pi4.or. + & tauangle(3,i).gt.pi34.and.tauangle(3,i).le.pi.or. + & tauangle(3,i).gt.-pi.and.tauangle(3,i).le.-pi34) then + call vecpr(dc_norm2(1,i-1+nres),dc_norm(1,i-2),vp1) + call vecpr(dc_norm2(1,i-2+nres),dc_norm2(1,i-1+nres),vp2) + call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3) + do j=1,3 + ctgt=cost/sint + ctgt1=cost1/sint1 + cosg_inv=1.0d0/cosg + dsintau(j,3,1,i)=-sing*ctgt1*domicron(j,2,1,i-1) + & -(fac0*vp1(j)+sing*dc_norm2(j,i-2+nres)) + & *vbld_inv(i-2+nres) + dtauangle(j,3,1,i)=cosg_inv*dsintau(j,3,1,i) + dsintau(j,3,2,i)= + & -sing*(ctgt1*domicron(j,2,2,i-1)+ctgt*domicron(j,1,1,i)) + & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1+nres) + dtauangle(j,3,2,i)=cosg_inv*dsintau(j,3,2,i) +c Bug fixed 3/24/05 (AL) + dsintau(j,3,3,i)=-sing*ctgt*domicron(j,1,2,i) + & +(fac0*vp3(j)-sing*dc_norm2(j,i-1+nres)) + & *vbld_inv(i-1+nres) +c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1) + dphi(j,3,3,i)=cosg_inv*dsintau(j,3,3,i) + enddo +c Obtaining the gamma derivatives from cosine derivative + else + do j=1,3 + dcostau(j,3,1,i)=fac1*dcosomicron(j,2,1,i-1)+fac3* + & dcostheta(j,1,i-1)-fac0*(dc_norm2(j,i-1+nres)-scalp* + & dc_norm2(j,i-2+nres))/vbld(i-2+nres) + dtauangle(j,3,1,i)=-1/sing*dcostau(j,3,1,i) + dcostau(j,3,2,i)=fac1*dcosomicron(j,2,2,i-1)+fac2* + & dcosomicron(j,1,1,i)+fac3*dcosomicron(j,2,2,i-1)+fac4* + & dcosomicron(j,1,1,i) + dtauangle(j,3,2,i)=-1/sing*dcostau(j,3,2,i) + dcostau(j,3,3,i)=fac2*dcosomicron(j,1,2,i)+fac4* + & dcostau(j,3,2,i)-fac0*(dc_norm2(j,i-2+nres)-scalp* + & dc_norm2(j,i-1+nres))/vbld(i-1+nres) + dtauangle(j,3,3,i)=-1/sing*dcostau(j,3,3,i) + enddo + endif + enddo #ifdef CRYST_SC c Derivatives of side-chain angles alpha and omega #if defined(MPI) && defined(PARINTDER) diff --git a/source/unres/src_MD/parmread.F b/source/unres/src_MD/parmread.F index be5f8b8..be1a193 100644 --- a/source/unres/src_MD/parmread.F +++ b/source/unres/src_MD/parmread.F @@ -553,33 +553,57 @@ C enddo endif #endif +C Read of Side-chain backbone correlation parameters +C Modified 11 May 2012 by Adasko +CCC C -C 5/21/07 (AL) Read coefficients of the backbone-local sidechain-local -C correlation energies. -C - read (isccor,*,end=119,err=119) nterm_sccor - do i=1,20 - do j=1,20 - read (isccor,'(a)') - do k=1,nterm_sccor - read (isccor,*,end=119,err=119) kk,v1sccor(k,i,j), - & v2sccor(k,i,j) + read (isccor,*,end=113,err=113) nsccortyp + read (isccor,*,end=113,err=113) (isccortyp(i),i=1,ntyp) +c write (iout,*) 'ntortyp',ntortyp + maxinter=3 +cc maxinter is maximum interaction sites + do l=1,maxinter + do i=1,nsccortyp + do j=1,nsccortyp + read (isccor,*,end=113,err=113) nterm_sccor(i,j),nlor_sccor(i,j) + v0ijsccor=0.0d0 + si=-1.0d0 + do k=1,nterm_sccor(i,j) + read (isccor,*,end=113,err=113) kk,v1sccor(k,i,j),v2sccor(k,i,j) + v0ijsccor=v0ijsccor+si*v1sccor(k,i,j) + si=-si + enddo + do k=1,nlor_sccor(i,j) + read (isccor,*,end=113,err=113) kk,vlor1sccor(k,i,j), + & vlor2sccor(k,i,j),vlor3sccor(k,i,j) + v0ijsccor=v0ijsccor+vlor1sccor(k,i,j)/ + &(1+vlor3sccor(k,i,j)**2) enddo + v0sccor(i,j)=v0ijsccor enddo enddo + enddo close (isccor) + if (lprint) then - write (iout,'(/a/)') 'Torsional constants of SCCORR:' - do i=1,20 - do j=1,20 + write (iout,'(/a/)') 'Torsional constants:' + do i=1,nsccortyp + do j=1,nsccortyp write (iout,*) 'ityp',i,' jtyp',j - do k=1,nterm_sccor + write (iout,*) 'Fourier constants' + do k=1,nterm_sccor(i,j) write (iout,'(2(1pe15.5))') v1sccor(k,i,j),v2sccor(k,i,j) enddo + write (iout,*) 'Lorenz constants' + do k=1,nlor_sccor(i,j) + write (iout,'(3(1pe15.5))') + & vlor1sccor(k,i,j),vlor2sccor(k,i,j),vlor3sccor(k,i,j) + enddo enddo enddo endif C +C C 9/18/99 (AL) Read coefficients of the Fourier expansion of the local C interaction energy of the Gly, Ala, and Pro prototypes. C -- 1.7.9.5