From: Adam Kazimierz Sieradzan Date: Thu, 14 Jun 2012 11:03:08 +0000 (-0400) Subject: Dzialajacy INTERTYP=1 i INTERTYP=2 X-Git-Tag: v.3.2~79^2~7 X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?p=unres.git;a=commitdiff_plain;h=ba3b4257f01c3857bbae7927ed7815ab7f5419d6 Dzialajacy INTERTYP=1 i INTERTYP=2 --- diff --git a/bin/unres/MD/unres_ifort_MPICH_GAB.exe b/bin/unres/MD/unres_ifort_MPICH_GAB.exe index f772e87..e333657 100755 Binary files a/bin/unres/MD/unres_ifort_MPICH_GAB.exe and b/bin/unres/MD/unres_ifort_MPICH_GAB.exe differ diff --git a/source/unres/src_MD/COMMON.CHAIN b/source/unres/src_MD/COMMON.CHAIN index f7a8a1d..6e19f8d 100644 --- a/source/unres/src_MD/COMMON.CHAIN +++ b/source/unres/src_MD/COMMON.CHAIN @@ -1,9 +1,10 @@ integer nres,nsup,nstart_sup,nz_start,nz_end,iz_sc, & nres0,nstart_seq double precision c,dc,dc_old,d_c_work,xloc,xrot,dc_norm,t,r, - & prod,rt,dc_work,cref,crefjlee + & prod,rt,dc_work,cref,crefjlee,dc_norm2 common /chain/ c(3,maxres2+2),dc(3,0:maxres2),dc_old(3,0:maxres2), & xloc(3,maxres),xrot(3,maxres),dc_norm(3,0:maxres2), + & dc_norm2(3,0:maxres2), & dc_work(MAXRES6),nres,nres0 common /rotmat/ t(3,3,maxres),r(3,3,maxres),prod(3,3,maxres), & rt(3,3,maxres) diff --git a/source/unres/src_MD/COMMON.SCCOR b/source/unres/src_MD/COMMON.SCCOR index 1c725d7..ca6210f 100644 --- a/source/unres/src_MD/COMMON.SCCOR +++ b/source/unres/src_MD/COMMON.SCCOR @@ -1,11 +1,16 @@ - Parameters of the SCCOR term +cc Parameters of the SCCOR term double precision v1sccor,v2sccor,vlor1sccor, - & vlor2sccor,vlor3sccor + & vlor2sccor,vlor3sccor,gloc_sc, + & dcostau,dsintau,dtauangle,dcosomicron, + & domicron integer nterm_sccor,isccortyp,nsccortyp,nlor_sccor - common/torsion/v1sccor(maxterm_sccor,20,20), - & v2sccor(maxterm_sccor,20,20), + common/sccor/v1sccor(maxterm_sccor,3,20,20), + & v2sccor(maxterm_sccor,3,20,20), + & v0sccor(maxterm_sccor,20), & 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) - + & vlor3sccor(maxterm_sccor,20,20),gloc_sc(3,0:maxres2,10), + & dcostau(3,3,3,maxres2),dsintau(3,3,3,maxres2), + & dtauangle(3,3,3,maxres2),dcosomicron(3,3,3,maxres2), + & domicron(3,3,3,maxres2) diff --git a/source/unres/src_MD/COMMON.VAR b/source/unres/src_MD/COMMON.VAR index c2a0a25..d560c87 100644 --- a/source/unres/src_MD/COMMON.VAR +++ b/source/unres/src_MD/COMMON.VAR @@ -5,7 +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) + & 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/DIMENSIONS b/source/unres/src_MD/DIMENSIONS index 224dade..5151ff7 100644 --- a/source/unres/src_MD/DIMENSIONS +++ b/source/unres/src_MD/DIMENSIONS @@ -54,7 +54,7 @@ C virtual-bond angle bending potentials & mmaxtheterm=maxtheterm) c Max number of torsional terms in SCCOR integer maxterm_sccor - parameter (maxterm_sccor=3) + parameter (maxterm_sccor=6) C Max. number of lobes in SC distribution integer maxlob parameter (maxlob=4) diff --git a/source/unres/src_MD/Makefile_MPICH_ifort b/source/unres/src_MD/Makefile_MPICH_ifort index e33762c..8846882 100644 --- a/source/unres/src_MD/Makefile_MPICH_ifort +++ b/source/unres/src_MD/Makefile_MPICH_ifort @@ -15,7 +15,7 @@ INSTALL_DIR = /users/software/mpich-1.2.7p1_intel-10.1_em64_ssh FC= ifort -OPT = -O3 -ip -w +OPT = -g -ip -w -CB FFLAGS = -c ${OPT} -I$(INSTALL_DIR)/include FFLAGS1 = -c -w -g -d2 -CA -CB -I$(INSTALL_DIR)/include diff --git a/source/unres/src_MD/checkder_p.F b/source/unres/src_MD/checkder_p.F index 719770e..4d0379e 100644 --- a/source/unres/src_MD/checkder_p.F +++ b/source/unres/src_MD/checkder_p.F @@ -8,6 +8,7 @@ C Check the gradient of Cartesian coordinates in internal coordinates. include 'COMMON.GEO' include 'COMMON.LOCAL' include 'COMMON.DERIV' + include 'COMMON.SCCOR' dimension temp(6,maxres),xx(3),gg(3) indmat(i,j)=((2*(nres-2)-i)*(i-1))/2+j-1 * @@ -180,6 +181,7 @@ C Check the gradient of the energy in Cartesian coordinates. include 'COMMON.IOUNITS' include 'COMMON.VAR' include 'COMMON.CONTACTS' + include 'COMMON.SCCOR' common /srutu/ icall dimension ggg(6),cc(3),xx(3),ddc(3),ddx(3),x(maxvar),g(maxvar) dimension grad_s(6,maxres) @@ -261,6 +263,7 @@ C Check the gradient of the energy in Cartesian coordinates. include 'COMMON.MD' include 'COMMON.LOCAL' include 'COMMON.SPLITELE' + include 'COMMON.SCCOR' common /srutu/ icall dimension ggg(6),ggg1(6),cc(3),xx(3),ddc(3),ddx(3),x(maxvar), & g(maxvar) @@ -288,10 +291,16 @@ c call checkintcartgrad call geom_to_var(nvar,x) if (.not.split_ene) then call etotal(energia(0)) +c do i=1,nres +c write (iout,*) "atu?", gloc_sc(1,i,icg),gloc(i,icg) +c enddo etot=energia(0) call enerprint(energia(0)) call flush(iout) write (iout,*) "enter cartgrad" +c do i=1,nres +c write (iout,*) gloc_sc(1,i,icg) +c enddo call flush(iout) call cartgrad write (iout,*) "exit cartgrad" @@ -338,6 +347,9 @@ c call checkintcartgrad call zerograd call etotal_short(energia(0)) call enerprint(energia(0)) +c do i=1,nres +c write (iout,*) gloc_sc(1,i,icg) +c enddo call flush(iout) write (iout,*) "enter cartgrad" call flush(iout) @@ -515,7 +527,7 @@ c------------------------------------------------------------------------- be=0.0D0 if (i.gt.2) then if (i.le.nres) phi(i+1)=beta(i-2,i-1,i,i+1) - if (itype(i).ne.10).and.(itype(i-1).ne.10) then + 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 diff --git a/source/unres/src_MD/cinfo.f b/source/unres/src_MD/cinfo.f index d3e68f7..cd1bd7b 100644 --- a/source/unres/src_MD/cinfo.f +++ b/source/unres/src_MD/cinfo.f @@ -1,11 +1,11 @@ C DO NOT EDIT THIS FILE - IT HAS BEEN GENERATED BY COMPINFO.C -C 2 5 38 +C 2 5 250 subroutine cinfo include 'COMMON.IOUNITS' write(iout,*)'++++ Compile info ++++' - write(iout,*)'Version 2.5 build 38' - write(iout,*)'compiled Mon Apr 2 08:28:34 2012' - write(iout,*)'compiled by adam@matrix.chem.cornell.edu' + write(iout,*)'Version 2.5 build 250' + write(iout,*)'compiled Thu Jun 14 07:01:36 2012' + write(iout,*)'compiled by aks255@matrix.chem.cornell.edu' write(iout,*)'OS name: Linux ' write(iout,*)'OS release: 2.6.34.9-69.fc13.x86_64 ' write(iout,*)'OS version:', @@ -13,7 +13,7 @@ C 2 5 38 write(iout,*)'flags:' write(iout,*)'INSTALL_DIR = /users/software/mpich-1.2.7p1_int...' write(iout,*)'FC= ifort' - write(iout,*)'OPT = -O3 -ip -w ' + write(iout,*)'OPT = -g -ip -w -CB ' write(iout,*)'FFLAGS = -c ${OPT} -I$(INSTALL_DIR)/include ' write(iout,*)'FFLAGS1 = -c -w -g -d2 -CA -CB -I$(INSTALL_DIR)...' write(iout,*)'FFLAGS2 = -c -w -g -O0 -I$(INSTALL_DIR)/include ' diff --git a/source/unres/src_MD/energy_p_new_barrier.F b/source/unres/src_MD/energy_p_new_barrier.F index 3257273..3979f64 100644 --- a/source/unres/src_MD/energy_p_new_barrier.F +++ b/source/unres/src_MD/energy_p_new_barrier.F @@ -483,6 +483,7 @@ cMS$ATTRIBUTES C :: proc_proc include 'COMMON.CONTROL' include 'COMMON.TIME1' include 'COMMON.MAXGRAD' + include 'COMMON.SCCOR' #ifdef TIMING #ifdef MPI time01=MPI_Wtime() @@ -755,7 +756,6 @@ c enddo & +wturn3*gel_loc_turn3(i) & +wturn6*gel_loc_turn6(i) & +wel_loc*gel_loc_loc(i) - & +wsccor*gsccor_loc(i) enddo #ifdef DEBUG write (iout,*) "gloc after adding corr" @@ -5650,7 +5650,7 @@ C Proline-Proline pair is a special case... & restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,itori,itori1, & (v1(j,itori,itori1),j=1,6),(v2(j,itori,itori1),j=1,6) gloc(i-3,icg)=gloc(i-3,icg)+wtor*gloci -c write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) + write (iout,*) 'i=',i,' gloc=',gloc(i-3,icg) enddo ! 6/20/98 - dihedral angle constraints edihcnstr=0.0d0 @@ -5765,6 +5765,7 @@ c do i=1,ndih_constr else difi=0.0 endif +c write (iout,*) "gloci", gloc(i-3,icg) cd write (iout,'(2i5,4f8.3,2e14.5)') i,itori,rad2deg*phii, cd & rad2deg*phi0(i), rad2deg*drange(i), cd & rad2deg*difi,0.25d0*ftors*difi**4,gloc(itori-3,icg) @@ -5835,6 +5836,7 @@ c lprn=.true. enddo gloc(i-3,icg)=gloc(i-3,icg)+wtor_d*gloci1 gloc(i-2,icg)=gloc(i-2,icg)+wtor_d*gloci2 +c write (iout,*) "gloci", gloc(i-3,icg) enddo return end @@ -5879,18 +5881,21 @@ cc Omicron is flat angle depending on the value of first digit c(see comment below) gloci=0.0D0 - do intertyp=1,3 + do intertyp=1,1 !intertyp 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 - if ((intertyp.eq.2).and.(i.eq.iphi_start-1)) cycle - if ((intertyp.eq.1).and.(i.eq.iphi_end+1)) cycle + if (((intertyp.eq.3).and.((itype(i-2).eq.10).or. + & (itype(i-1).eq.10).or.(itype(i-2).eq.21).or. + & (itype(i-1).eq.21))) + & .or. ((intertyp.eq.1).and.((itype(i-2).eq.10) + & .or.(itype(i-2).eq.21))) + & .or.((intertyp.eq.2).and.((itype(i-1).eq.10).or. + & (itype(i-1).eq.21)))) cycle + if ((intertyp.eq.2).and.(i.le.iphi_start-1)) cycle + if ((intertyp.eq.1).and.(i.ge.iphi_end+1)) cycle do j=1,nterm_sccor(isccori,isccori1) v1ij=v1sccor(j,intertyp,isccori,isccori1) v2ij=v2sccor(j,intertyp,isccori,isccori1) @@ -5900,14 +5905,19 @@ c 3 = SC...Ca...Ca...SC gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi) enddo gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci +c write (iout,*) "WTF",intertyp,i,itype(i), +c &gloc_sc(intertyp,i-3,icg) 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,intertyp,itori,itori1),j=1,6) & ,(v2sccor(j,intertyp,itori,itori1),j=1,6) gsccor_loc(i-3)=gsccor_loc(i-3)+gloci + enddo !intertyp enddo - enddo +c do i=1,nres +c write (iout,*) "W@T@F", gloc_sc(1,i,icg),gloc(i,icg) +c enddo return end c---------------------------------------------------------------------------- diff --git a/source/unres/src_MD/gradient_p.F b/source/unres/src_MD/gradient_p.F index 375fcf4..7fec1e8 100644 --- a/source/unres/src_MD/gradient_p.F +++ b/source/unres/src_MD/gradient_p.F @@ -8,6 +8,7 @@ include 'COMMON.FFIELD' include 'COMMON.MD' include 'COMMON.IOUNITS' + include 'COMMON.SCCOR' external ufparm integer uiparm(1) double precision urparm(1) @@ -263,10 +264,15 @@ C------------------------------------------------------------------------- include 'COMMON.MD' include 'COMMON.IOUNITS' include 'COMMON.TIME1' + include 'COMMON.SCCOR' c c This subrouting calculates total Cartesian coordinate gradient. c The subroutine chainbuild_cart and energy MUST be called beforehand. c +c do i=1,nres +c write (iout,*) "przed sum_grad", gloc_sc(1,i,icg),gloc(i,icg) +c enddo + #ifdef TIMING time00=MPI_Wtime() #endif @@ -274,6 +280,9 @@ c call sum_gradient #ifdef TIMING #endif +c do i=1,nres +c write (iout,*) "checkgrad", gloc_sc(1,i,icg),gloc(i,icg) +c enddo cd write (iout,*) "After sum_gradient" cd do i=1,nres-1 cd write (iout,*) i," gradc ",(gradc(j,i,icg),j=1,3) @@ -337,6 +346,7 @@ C------------------------------------------------------------------------- include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.MD' + include 'COMMON.SCCOR' C C Initialize Cartesian-coordinate gradient C @@ -374,6 +384,9 @@ C gradx(j,i,icg)=0.0d0 gscloc(j,i)=0.0d0 gsclocx(j,i)=0.0d0 + do intertyp=1,3 + gloc_sc(intertyp,i,icg)=0.0d0 + enddo enddo enddo C diff --git a/source/unres/src_MD/int_to_cart.f b/source/unres/src_MD/int_to_cart.f index 149e86f..b2d378e 100644 --- a/source/unres/src_MD/int_to_cart.f +++ b/source/unres/src_MD/int_to_cart.f @@ -13,9 +13,15 @@ c------------------------------------------------------------- include 'COMMON.INTERACT' include 'COMMON.MD' include 'COMMON.IOUNITS' - + include 'COMMON.SCCOR' c calculating dE/ddc1 if (nres.lt.3) goto 18 +c do i=1,nres +c c do intertyp=1,3 +c write (iout,*) "przed tosyjnymi",i,intertyp,gcart(intertyp,i) +c &,gloc_sc(1,i,icg),gloc(i,icg) +c enddo +c enddo do j=1,3 gcart(j,1)=gcart(j,1)+gloc(1,icg)*dphi(j,1,4) & +gloc(nres-2,icg)*dtheta(j,1,3) @@ -119,9 +125,13 @@ C INTERTYP=2 Ca...Ca...Ca...SC C INTERTYP=3 SC...Ca...Ca...SC c calculating dE/ddc1 18 continue +c do i=1,nres +c gloc(i,icg)=0.0D0 +c write (iout,*) "poczotkoawy",i,gloc_sc(1,i,icg) +c enddo if (nres.lt.2) return if ((nres.lt.3).and.(itype(1).eq.10)) return - if (itype(1).ne.10) then + if ((itype(1).ne.10).and.(itype(1).ne.21)) then do j=1,3 cc Derviative was calculated for oposite vector of side chain therefore c there is "-" sign before gloc_sc @@ -129,15 +139,16 @@ c there is "-" sign before gloc_sc & dtauangle(j,1,1,3) gcart(j,1)=gcart(j,1)+gloc_sc(1,0,icg)* & dtauangle(j,1,2,3) - if (itype(2).ne.10) then + if ((itype(2).ne.10).and.(itype(2).ne.21)) then gxcart(j,1)= gxcart(j,1) & -gloc_sc(3,0,icg)*dtauangle(j,3,1,3) gcart(j,1)=gcart(j,1)+gloc_sc(3,0,icg)* - dtauangle(j,3,2,3) + & dtauangle(j,3,2,3) endif enddo endif - if (nres.ge.3).and.(itype(3).ne.10) then + if ((nres.ge.3).and.(itype(3).ne.10).and.(itype(3).ne.21)) + & then do j=1,3 gcart(j,1)=gcart(j,1)+gloc_sc(2,1,icg)*dtauangle(j,2,1,4) enddo @@ -148,36 +159,51 @@ c & +gloc_sc(intertyp,nres-2,icg)*dtheta(j,1,3) c Calculating the remainder of dE/ddc2 do j=1,3 - if(itype(2).ne.10) then - if (itype(1).ne.10) gxcart=(j,2)=gxcart(j,2)+ + if((itype(2).ne.10).and.(itype(2).ne.21)) then + if (itype(1).ne.10) gxcart(j,2)=gxcart(j,2)+ & gloc_sc(3,0,icg)*dtauangle(j,3,3,3) - if ((itype(3).ne.10).and.(nres.ge.3)) then + if ((itype(3).ne.10).and.(nres.ge.3).and.(itype(3).ne.21)) then gxcart(j,2)=gxcart(j,2)-gloc_sc(3,1,icg)*dtauangle(j,3,1,4) cc the - above is due to different vector direction gcart(j,2)=gcart(j,2)+gloc_sc(3,1,icg)*dtauangle(j,3,2,4) endif if (nres.gt.3) then - gxcart(j,2)=gxcart(j,2)-gloc_sc(2,1,icg)*dtauangle(j,1,1,4) + gxcart(j,2)=gxcart(j,2)-gloc_sc(1,1,icg)*dtauangle(j,1,1,4) cc the - above is due to different vector direction - gcart(j,2)=gcart(j,2)+gloc_sc(2,1,icg)*dtauangle(j,1,2,4) + gcart(j,2)=gcart(j,2)+gloc_sc(1,1,icg)*dtauangle(j,1,2,4) +c write(iout,*) gloc_sc(1,1,icg),dtauangle(j,1,2,4),"gcart" +c write(iout,*) gloc_sc(1,1,icg),dtauangle(j,1,1,4),"gx" endif endif - if (itype(1).ne.10) then + if ((itype(1).ne.10).and.(itype(1).ne.21)) then gcart(j,2)=gcart(j,2)+gloc_sc(1,0,icg)*dtauangle(j,1,3,3) +c write(iout,*) gloc_sc(1,0,icg),dtauangle(j,1,3,3) endif - if ((itype(3).ne.10).and(nres.ge.3)) then + if ((itype(3).ne.10).and.(nres.ge.3)) then gcart(j,2)=gcart(j,2)+gloc_sc(2,1,icg)*dtauangle(j,2,2,4) + write(iout,*) gloc_sc(2,1,icg),dtauangle(j,2,2,4) + endif + if ((itype(4).ne.10).and.(nres.ge.4)) then + gcart(j,2)=gcart(j,2)+gloc_sc(2,2,icg)*dtauangle(j,2,1,5) +c write(iout,*) gloc_sc(2,2,icg),dtauangle(j,2,1,5) endif + +c write(iout,*) gcart(j,2),itype(2),itype(1),itype(3), "gcart2" + enddo c If there are more than five residues if(nres.ge.5) then do i=3,nres-2 do j=1,3 +c write(iout,*) "before", gcart(j,i) if (itype(i).ne.10) then gxcart(j,i)=gxcart(j,i)+gloc_sc(2,i-2,icg) & *dtauangle(j,2,3,i+1) - & -gloc_sc(1,i-1,icg)*dtauangle(j,1,1i+2) - gcart(j,i)=gcart(j,i)+gloc_sc(1,i-1,icg)* + & -gloc_sc(1,i-1,icg)*dtauangle(j,1,1,i+2) + gcart(j,i)=gcart(j,i)+gloc_sc(1,i-1,icg) & *dtauangle(j,1,2,i+2) +c write(iout,*) "new",j,i, +c & gcart(j,i),gloc_sc(1,i-1,icg),dtauangle(j,1,2,i+2) + if (itype(i-1).ne.10) then gxcart(j,i)=gxcart(j,i)+gloc_sc(3,i-2,icg) &*dtauangle(j,3,3,i+1) @@ -196,6 +222,8 @@ c If there are more than five residues if (itype(i+1).ne.10) then gcart(j,i)=gcart(j,i)+gloc_sc(2,i-1,icg)* & dtauangle(j,2,2,i+2) +c write(iout,*) "numer",i,gloc_sc(2,i-1,icg), +c & dtauangle(j,2,2,i+2) endif if (itype(i+2).ne.10) then gcart(j,i)=gcart(j,i)+gloc_sc(2,i,icg)* @@ -207,37 +235,42 @@ c If there are more than five residues c Setting dE/ddnres-1 if(nres.ge.4) then do j=1,3 - if (itype(nres-1).ne.10) then + if ((itype(nres-1).ne.10).and.(itype(nres-1).ne.21)) then gxcart(j,nres-1)=gxcart(j,nres-1)+gloc_sc(2,nres-3,icg) & *dtauangle(j,2,3,nres) + write (iout,*) "gxcart(nres-1)", gloc_sc(2,nres-3,icg), + & dtauangle(j,2,3,nres), gxcart(j,nres-1) if (itype(nres-2).ne.10) then gxcart(j,nres-1)=gxcart(j,nres-1)+gloc_sc(3,nres-3,icg) & *dtauangle(j,3,3,nres) endif - if (itype(nres).ne.10) then + if ((itype(nres).ne.10).and.(itype(nres).ne.21)) then gxcart(j,nres-1)=gxcart(j,nres-1)-gloc_sc(3,nres-2,icg) & *dtauangle(j,3,1,nres+1) gcart(j,nres-1)=gcart(j,nres-1)+gloc_sc(3,nres-2,icg) & *dtauangle(j,3,2,nres+1) endif endif - if (itype(nres-2).ne.10) then - gcart(j,nres-1)=gcart(j,nres-1)+gloc_sc(2,nres-3,icg)* - & *dtauangle(j,2,3,nres) + if ((itype(nres-2).ne.10).and.(itype(nres-2).ne.21)) then + gcart(j,nres-1)=gcart(j,nres-1)+gloc_sc(1,nres-3,icg)* + & dtauangle(j,1,3,nres) endif - if (itype(nres).ne.10) then + if ((itype(nres).ne.10).and.(itype(nres).ne.21)) then gcart(j,nres-1)=gcart(j,nres-1)+gloc_sc(2,nres-2,icg)* - *dtauangle(j,2,2,nres+1) + & dtauangle(j,2,2,nres+1) + write (iout,*) "gcart(nres-1)", gloc_sc(2,nres-2,icg), + & dtauangle(j,2,2,nres+1), itype(nres-1),itype(nres) endif enddo endif c Settind dE/ddnres - if (nres.ge.3).and(itype(nres).ne.10)then + if ((nres.ge.3).and.(itype(nres).ne.10))then do j=1,3 gxcart(j,nres)=gxcart(j,nres)+gloc_sc(3,nres-2,icg) & *dtauangle(j,3,3,nres+1)+gloc_sc(2,nres-2,icg) & *dtauangle(j,2,3,nres+1) enddo + endif c The side-chain vector derivatives return end diff --git a/source/unres/src_MD/intcartderiv.F b/source/unres/src_MD/intcartderiv.F index 8f4a71e..2068031 100644 --- a/source/unres/src_MD/intcartderiv.F +++ b/source/unres/src_MD/intcartderiv.F @@ -12,6 +12,7 @@ include 'COMMON.DERIV' include 'COMMON.IOUNITS' include 'COMMON.LOCAL' + include 'COMMON.SCCOR' double precision dcostheta(3,2,maxres), & dcosphi(3,3,maxres),dsinphi(3,3,maxres), & dcosalpha(3,3,maxres),dcosomega(3,3,maxres), @@ -51,20 +52,20 @@ c We need dtheta(:,:,i-1) to compute dphi(:,:,i) #else do i=3,nres #endif - if (itype(i-1).ne.10) then + if ((itype(i-1).ne.10).and.(itype(i-1).ne.21)) 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)+ + 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) + & +cost1*(dc_norm(j,i-1+nres)))/ + & vbld(i-1+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 @@ -73,8 +74,9 @@ CC Looks messy but better than if in loop & 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) + & +cost2*(-dc_norm(j,i-1+nres)))/ + & vbld(i-1+nres) +c write(iout,*) "vbld", i,itype(i),vbld(i-1+nres) domicron(j,2,2,i)=-1/sint2*dcosomicron(j,2,2,i) enddo endif @@ -153,6 +155,7 @@ Calculate derivative of Tauangle #else do i=3,nres #endif + if ((itype(i-2).eq.21).or.(itype(i-2).eq.10)) cycle cc dtauangle(j,intertyp,dervityp,residue number) cc INTERTYP=1 SC...Ca...Ca..Ca c the conventional case @@ -164,6 +167,7 @@ c the conventional case cosg=dcos(tauangle(1,i)) do j=1,3 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres) +cc write(iout,*) dc_norm2(j,i-2+nres),"dcnorm" enddo scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1)) fac0=1.0d0/(sint1*sint) @@ -171,13 +175,11 @@ c the conventional case fac2=cost1*fac0 fac3=cosg*cost1/(sint1*sint1) fac4=cosg*cost/(sint*sint) +cc write(iout,*) "faki",fac0,fac1,fac2,fac3,fac4 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) @@ -185,37 +187,42 @@ c Obtaining the gamma derivatives from sine derivative 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))) + dsintau(j,1,1,i)=-sing*ctgt1*domicron(j,2,2,i-1) + &-(fac0*vp1(j)+sing*(dc_norm2(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)) + & -sing*(ctgt1*domicron(j,2,1,i-1)+ctgt*dtheta(j,1,i)) & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1) +c write(iout,*) "dsintau", dsintau(j,1,2,i) 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) + dtauangle(j,1,3,i)=cosg_inv*dsintau(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) + dcostau(j,1,1,i)=fac1*dcosomicron(j,2,2,i-1)+fac3* + & dcosomicron(j,2,2,i-1)-fac0*(dc_norm(j,i-1)-scalp* + & (dc_norm2(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* + dcostau(j,1,2,i)=fac1*dcosomicron(j,2,1,i-1)+fac2* + & dcostheta(j,1,i)+fac3*dcosomicron(j,2,1,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) +c write (iout,*) "else",i enddo - endif + endif +c do k=1,3 +c write(iout,*) "tu",i,k,(dtauangle(j,1,k,i),j=1,3) +c enddo enddo CC Second case Ca...Ca...Ca...SC #ifdef PARINTDER @@ -223,11 +230,12 @@ CC Second case Ca...Ca...Ca...SC #else do i=4,nres #endif + if ((itype(i-1).eq.21).or.(itype(i-1).eq.10)) cycle c the conventional case sint=dsin(omicron(1,i)) sint1=dsin(theta(i-1)) - sing=dsin(tauangle(2,i) - cost=dcos(omicron(1,i) + sing=dsin(tauangle(2,i)) + cost=dcos(omicron(1,i)) cost1=dcos(theta(i-1)) cosg=dcos(tauangle(2,i)) c do j=1,3 @@ -251,17 +259,22 @@ c Obtaining the gamma derivatives from sine derivative 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) + & +(fac0*vp1(j)-sing*dc_norm(j,i-3))*vbld_inv(i-2) +c write(iout,*) i,j,dsintau(j,2,1,i),sing*ctgt1*dtheta(j,1,i-1), +c &fac0*vp1(j),sing*dc_norm(j,i-3),vbld_inv(i-2),"dsintau(2,1)" 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 write(iout,*) "sprawdzenie",i,j,sing*ctgt1*dtheta(j,2,i-1), +c & sing*ctgt*domicron(j,1,2,i), +c & (fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1) + dtauangle(j,2,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_norm(j,i-1+nres))*vbld_inv(i) + & +(fac0*vp3(j)-sing*dc_norm(j,i-1+nres))*vbld_inv(i-1+nres) 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) + dtauangle(j,2,3,i)=cosg_inv*dsintau(j,2,3,i) enddo c Obtaining the gamma derivatives from cosine derivative else @@ -273,11 +286,12 @@ c Obtaining the gamma derivatives from cosine derivative 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) + dtauangle(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* + & dcosomicron(j,1,2,i)-fac0*(dc_norm(j,i-3)-scalp* & dc_norm(j,i-1+nres))/vbld(i-1+nres) - dtauanlge(j,2,3,i)=-1/sing*dcosphi(j,3,i) + dtauangle(j,2,3,i)=-1/sing*dcostau(j,2,3,i) + write(iout,*) i,j,"else", dtauangle(j,2,3,i) enddo endif enddo @@ -291,6 +305,8 @@ CCC third case SC...Ca...Ca...SC do i=3,nres #endif c the conventional case + if ((itype(i-1).eq.21).or.(itype(i-1).eq.10).or. + &(itype(i-2).eq.21).or.(itype(i-2).eq.10)) cycle sint=dsin(omicron(1,i)) sint1=dsin(omicron(2,i-1)) sing=dsin(tauangle(3,i)) @@ -331,7 +347,7 @@ c Bug fixed 3/24/05 (AL) & +(fac0*vp3(j)-sing*dc_norm(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) + dtauangle(j,3,3,i)=cosg_inv*dsintau(j,3,3,i) enddo c Obtaining the gamma derivatives from cosine derivative else diff --git a/source/unres/src_MD/parmread.F b/source/unres/src_MD/parmread.F index be1a193..1968235 100644 --- a/source/unres/src_MD/parmread.F +++ b/source/unres/src_MD/parmread.F @@ -562,15 +562,17 @@ C c write (iout,*) 'ntortyp',ntortyp maxinter=3 cc maxinter is maximum interaction sites - do l=1,maxinter + 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) + read (isccor,*,end=113,err=113) kk,v1sccor(k,l,i,j) + & ,v2sccor(k,l,i,j) + v0ijsccor=v0ijsccor+si*v1sccor(k,l,i,j) si=-si enddo do k=1,nlor_sccor(i,j) @@ -592,7 +594,7 @@ cc maxinter is maximum interaction sites write (iout,*) 'ityp',i,' jtyp',j write (iout,*) 'Fourier constants' do k=1,nterm_sccor(i,j) - write (iout,'(2(1pe15.5))') v1sccor(k,i,j),v2sccor(k,i,j) + write (iout,'(2(1pe15.5))') v1sccor(k,l,i,j),v2sccor(k,l,i,j) enddo write (iout,*) 'Lorenz constants' do k=1,nlor_sccor(i,j) diff --git a/source/unres/src_MD/tmptmp b/source/unres/src_MD/tmptmp index 54e7a36..fb6f043 100644 --- a/source/unres/src_MD/tmptmp +++ b/source/unres/src_MD/tmptmp @@ -1 +1 @@ -adam +aks255