integer nres,nsup,nstart_sup,nz_start,nz_end,iz_sc,
& nres0,nstart_seq,chain_length,iprzes,tabperm,nperm
double precision c,dc,dc_old,d_c_work,xloc,xrot,dc_norm,t,r,
- & prod,rt,dc_work,cref,crefjlee,chain_rep
+ & prod,rt,dc_work,cref,crefjlee,chain_rep,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)
& iphi_end,iphid_start,iphid_end,ibond_start,ibond_end,
& ibondp_start,ibondp_end,ivec_start,ivec_end,iset_start,iset_end,
& iturn3_start,iturn3_end,iturn4_start,iturn4_end,iint_start,
- & iint_end,iphi1_start,iphi1_end,
+ & iint_end,iphi1_start,iphi1_end,itau_start,itau_end,
& ibond_displ(0:max_fg_procs-1),ibond_count(0:max_fg_procs-1),
& ithet_displ(0:max_fg_procs-1),ithet_count(0:max_fg_procs-1),
& iphi_displ(0:max_fg_procs-1),iphi_count(0:max_fg_procs-1),
& ibondp_start,ibondp_end,ivec_start,ivec_end,iset_start,iset_end,
& iturn3_start,iturn3_end,iturn4_start,iturn4_end,iint_start,
& iint_end,iphi1_start,iphi1_end,iint_count,iint_displ,ivec_displ,
- & ivec_count,iset_displ,
+ & ivec_count,iset_displ,itau_start,itau_end,
& iset_count,ibond_displ,ibond_count,ithet_displ,ithet_count,
& iphi_displ,iphi_count,iphi1_displ,iphi1_count
C Inverses of the actual virtual bond lengths
-C Parameters of the SCCOR term
- double precision v1sccor,v2sccor
- integer nterm_sccor
- common/sccor/v1sccor(maxterm_sccor,20,20),
- & v2sccor(maxterm_sccor,20,20),
- & nterm_sccor
+cc Parameters of the SCCOR term
+ double precision v1sccor,v2sccor,vlor1sccor,
+ & vlor2sccor,vlor3sccor,gloc_sc,
+ & dcostau,dsintau,dtauangle,dcosomicron,
+ & domicron
+ integer nterm_sccor,isccortyp,nsccortyp,nlor_sccor
+ common/sccor/v1sccor(maxterm_sccor,3,-ntyp:ntyp,-ntyp:ntyp),
+ & v2sccor(maxterm_sccor,3,-ntyp:ntyp,-ntyp:ntyp),
+ & v0sccor(maxterm_sccor,-ntyp:ntyp),
+ & nterm_sccor(-ntyp:ntyp,-ntyp:ntyp),isccortyp(-ntyp:ntyp),
+ & nsccortyp,
+ & nlor_sccor(-ntyp:ntyp,-ntyp:ntyp),
+ & vlor1sccor(maxterm_sccor,20,20),
+ & vlor2sccor(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)
& mask_theta,mask_phi,mask_side
double precision theta,phi,alph,omeg,varsave,esave,varall,vbld,
& thetaref,phiref,costtab,sinttab,cost2tab,sint2tab,
+ & tauangle,omicron,
& xxtab,yytab,zztab,xxref,yyref,zzref
common /var/ theta(maxres),phi(maxres),alph(maxres),omeg(maxres),
& vbld(2*maxres),thetaref(maxres),phiref(maxres),
& costtab(maxres), sinttab(maxres), cost2tab(maxres),
+ & omicron(2,maxres),tauangle(3,maxres),
& sint2tab(maxres),xxtab(maxres),yytab(maxres),
& zztab(maxres),xxref(maxres),yyref(maxres),zzref(maxres),
& ialph(maxres,2),ivar(4*maxres2),ntheta,nphi,nside,nvar
& 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)
stdforcp(i)=stdfp*dsqrt(gamp)
enddo
do i=nnt,nct
- stdforcsc(i)=stdfsc(itype(i))*dsqrt(gamsc(itype(i)))
+ stdforcsc(i)=stdfsc(itype(i))*dsqrt(gamsc(iabs(itype(i))))
enddo
endif
c Open the pdb file for snapshotshots
+++ /dev/null
-Makefile_MPICH_ifort
\ No newline at end of file
--- /dev/null
+CPPFLAGS = -DLINUX -DUNRES -DMP -DMPI \
+ -DPGI -DSPLITELE -DISNAN -DAMD64 \
+ -DPROCOR -DLANG0 \
+ -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORPDB
+## -DPROCOR
+## -DMOMENT
+#-DCO_BIAS
+#-DCRYST_TOR
+#-DDEBUG
+
+#INSTALL_DIR = /usr/local/mpich-1.2.0
+INSTALL_DIR = /users/software/mpich-1.2.7p1_intel-10.1_em64_ssh
+#
+#FC= /usr/local/opt/intel/compiler60/ia32/bin/ifc
+FC= ifort
+
+#OPT = -O3 -ip -w
+OPT = -g -CB
+#OPT = -g
+CFLAGS = -DSGI -c
+
+FFLAGS = -c ${OPT} -I$(INSTALL_DIR)/include
+FFLAGS1 = -c -w -g -d2 -CA -CB -I$(INSTALL_DIR)/include
+FFLAGS2 = -c -w -O0 -I$(INSTALL_DIR)/include
+FFLAGSE = -c -w -O3 -ipo -ipo_obj -opt_report -I$(INSTALL_DIR)/include
+
+BIN = ../../../bin/unres/MD-M/unres_Tc_procor_oldparm_em64-D-symetr.exe
+#LIBS = -L$(INSTALL_DIR)/lib_pgi -lmpich xdrf/libxdrf.a
+#LIBS = -L$(INSTALL_DIR)/lib_ifort -lmpich xdrf/libxdrf.a
+LIBS = -L$(INSTALL_DIR)/lib -lmpich ../../lib/xdrf_em64/libxdrf.a -g -d2 -CA -CB
+
+ARCH = LINUX
+PP = /lib/cpp -P
+
+
+all: unres
+
+.SUFFIXES: .F
+.F.o:
+ ${FC} ${FFLAGS} ${CPPFLAGS} $*.F
+
+
+object = unres.o arcos.o cartprint.o chainbuild.o convert.o initialize_p.o \
+ matmult.o readrtns_CSA.o parmread.o gen_rand_conf.o printmat.o map.o \
+ pinorm.o randgens.o rescode.o intcor.o timing.o misc.o intlocal.o \
+ cartder.o checkder_p.o econstr_local.o energy_p_new_barrier.o \
+ energy_p_new-sep_barrier.o gradient_p.o minimize_p.o sumsld.o \
+ cored.o rmdd.o geomout.o readpdb.o permut.o regularize.o thread.o fitsq.o mcm.o \
+ mc.o bond_move.o refsys.o check_sc_distr.o check_bond.o contact.o djacob.o \
+ eigen.o blas.o add.o entmcm.o minim_mcmf.o \
+ together.o csa.o minim_jlee.o shift.o diff12.o bank.o newconf.o ran.o \
+ indexx.o MP.o compare_s1.o prng_32.o \
+ test.o banach.o distfit.o rmsd.o elecont.o dihed_cons.o \
+ sc_move.o local_move.o \
+ intcartderiv.o lagrangian_lesyng.o\
+ stochfric.o kinetic_lesyng.o MD_A-MTS.o moments.o int_to_cart.o \
+ surfatom.o sort.o muca_md.o MREMD.o rattle.o gauss.o energy_split-sep.o \
+ q_measure.o gnmr1.o
+
+unres: ${object} proc_proc.o
+ cc -o compinfo compinfo.c
+ ./compinfo | true
+ ${FC} ${FFLAGS} cinfo.f
+ ${FC} ${OPT} -Wl,-M ${object} proc_proc.o cinfo.o ${LIBS} -o ${BIN}
+
+
+clean:
+ /bin/rm *.o
+
+newconf.o: newconf.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} newconf.f
+
+bank.o: bank.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} bank.F
+
+diff12.o: diff12.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} diff12.f
+
+csa.o: csa.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} csa.f
+
+shift.o: shift.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} shift.F
+
+ran.o: ran.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} ran.f
+
+together.o: together.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} together.F
+
+test.o: test.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} test.F
+
+chainbuild.o: chainbuild.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} chainbuild.F
+
+matmult.o: matmult.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} matmult.f
+
+parmread.o : parmread.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} parmread.F
+
+intcor.o : intcor.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} intcor.f
+
+cartder.o : cartder.F
+ ${FC} ${FFLAGS} ${CPPFLAGS} cartder.F
+
+readpdb.o : readpdb.F
+ ${FC} ${FFLAGS2} ${CPPFLAGS} readpdb.F
+
+permut.o : permut.F
+ ${FC} ${FFLAGS2} ${CPPFLAGS} permut.F
+
+sumsld.o : sumsld.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} sumsld.f
+
+cored.o : cored.f
+ ${FC} ${FFLAGS1} ${CPPFLAGS} cored.f
+
+rmdd.o : rmdd.f
+ ${FC} ${FFLAGS} ${CPPFLAGS} rmdd.f
+
+energy_p_new.o : energy_p_new.F
+ ${FC} ${FFLAGSE} ${CPPFLAGS} energy_p_new.F
+
+lagrangian_lesyng.o : lagrangian_lesyng.F
+ ${FC} ${FFLAGSE} ${CPPFLAGS} lagrangian_lesyng.F
+
+proc_proc.o: proc_proc.c
+ ${CC} ${CFLAGS} proc_proc.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
+ 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
+ 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)
--- /dev/null
+C DO NOT EDIT THIS FILE - IT HAS BEEN GENERATED BY COMPINFO.C
+C 2 3 3438
+ subroutine cinfo
+ include 'COMMON.IOUNITS'
+ write(iout,*)'++++ Compile info ++++'
+ write(iout,*)'Version 2.3 build 3438'
+ write(iout,*)'compiled Thu Oct 4 09:08:03 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:',
+ & ' #1 SMP Tue May 3 09:23:03 UTC 2011 '
+ write(iout,*)'flags:'
+ write(iout,*)'CPPFLAGS = -DLINUX -DUNRES -DMP -DMPI \\'
+ write(iout,*)' -DPGI -DSPLITELE -DISNAN -DAMD64 \\'
+ write(iout,*)' -DPROCOR -DLANG0 \\'
+ write(iout,*)' -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORP...'
+ write(iout,*)'INSTALL_DIR = /users/software/mpich-1.2.7p1_int...'
+ write(iout,*)'FC= ifort'
+ write(iout,*)'OPT = -g -CB'
+ write(iout,*)'CFLAGS = -DSGI -c'
+ 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 -O0 -I$(INSTALL_DIR)/include'
+ write(iout,*)'FFLAGSE = -c -w -O3 -ipo -ipo_obj -opt_report ...'
+ write(iout,*)'BIN = ../../../bin/unres/MD-M/unres_Tc_procor_o...'
+ write(iout,*)'LIBS = -L$(INSTALL_DIR)/lib -lmpich ../../lib/x...'
+ write(iout,*)'ARCH = LINUX'
+ write(iout,*)'PP = /lib/cpp -P'
+ write(iout,*)'object = unres.o arcos.o cartprint.o chainbuild...'
+ write(iout,*)'++++ End of compile info ++++'
+ return
+ end
include 'mpif.h'
#endif
double precision gradbufc(3,maxres),gradbufx(3,maxres),
- & glocbuf(4*maxres),gradbufc_sum(3,maxres)
+ & glocbuf(4*maxres),gradbufc_sum(3,maxres),gloc_scbuf(3,maxres)
+#endif
include 'COMMON.SETUP'
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.CONTROL'
include 'COMMON.TIME1'
include 'COMMON.MAXGRAD'
+ include 'COMMON.SCCOR'
#ifdef TIMING
time01=MPI_Wtime()
#endif
& +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"
do i=1,4*nres
glocbuf(i)=gloc(i,icg)
enddo
+#define DEBUG
+#ifdef DEBUG
+ write (iout,*) "gloc_sc before reduce"
+ do i=1,nres
+ do j=1,1
+ write (iout,*) i,j,gloc_sc(j,i,icg)
+ enddo
+ enddo
+#endif
+#undef DEBUG
+ do i=1,nres
+ do j=1,3
+ gloc_scbuf(j,i)=gloc_sc(j,i,icg)
+ enddo
+ enddo
time00=MPI_Wtime()
call MPI_Barrier(FG_COMM,IERR)
time_barrier_g=time_barrier_g+MPI_Wtime()-time00
call MPI_Reduce(glocbuf(1),gloc(1,icg),4*nres,
& MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
time_reduce=time_reduce+MPI_Wtime()-time00
+ call MPI_Reduce(gloc_scbuf(1,1),gloc_sc(1,1,icg),3*nres,
+ & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
+ time_reduce=time_reduce+MPI_Wtime()-time00
+#define DEBUG
+#ifdef DEBUG
+ write (iout,*) "gloc_sc after reduce"
+ do i=1,nres
+ do j=1,1
+ write (iout,*) i,j,gloc_sc(j,i,icg)
+ enddo
+ enddo
+#endif
+#undef DEBUG
#ifdef DEBUG
write (iout,*) "gloc after reduce"
do i=1,4*nres
C Set lprn=.true. for debugging
lprn=.false.
c lprn=.true.
-c write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor
+c write (iout,*) "EBACK_SC_COR",itau_start,itau_end
esccor=0.0D0
- do i=iphi_start,iphi_end
- if (itype(i-2).eq.21 .or. itype(i-1).eq.21) cycle
+ do i=itau_start,itau_end
esccor_ii=0.0D0
- itori=itype(i-2)
- itori1=itype(i-1)
+ isccori=isccortyp(itype(i-2))
+ isccori1=isccortyp(itype(i-1))
+c write (iout,*) "EBACK_SC_COR",i,nterm_sccor(isccori,isccori1)
phii=phi(i)
+ do intertyp=1,3 !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...SCi
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)
+ if (((intertyp.eq.3).and.((itype(i-2).eq.10).or.
+ & (itype(i-1).eq.10).or.(itype(i-2).eq.ntyp1).or.
+ & (itype(i-1).eq.ntyp1)))
+ & .or. ((intertyp.eq.1).and.((itype(i-2).eq.10)
+ & .or.(itype(i-2).eq.ntyp1).or.(itype(i-1).eq.ntyp1)
+ & .or.(itype(i).eq.ntyp1)))
+ & .or.((intertyp.eq.2).and.((itype(i-1).eq.10).or.
+ & (itype(i-1).eq.ntyp1).or.(itype(i-2).eq.ntyp1).or.
+ & (itype(i-3).eq.ntyp1)))) cycle
+ if ((intertyp.eq.2).and.(i.eq.4).and.(itype(1).eq.ntyp1)) cycle
+ if ((intertyp.eq.1).and.(i.eq.nres).and.(itype(nres).eq.ntyp1))
+ & 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
+c write (iout,*) "EBACK_SC_COR",i,esccor,intertyp
+ gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*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)
+ & restyp(itype(i-2)),i-2,restyp(itype(i-1)),i-1,isccori,isccori1,
+ & (v1sccor(j,intertyp,isccori,isccori1),j=1,6)
+ & ,(v2sccor(j,intertyp,isccori,isccori1),j=1,6)
gsccor_loc(i-3)=gsccor_loc(i-3)+gloci
+ enddo !intertyp
enddo
+
return
end
c----------------------------------------------------------------------------
include 'COMMON.CHAIN'
include 'COMMON.VAR'
include 'COMMON.MD'
+ include 'COMMON.SCCOR'
C
C Initialize Cartesian-coordinate gradient
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
iphi_end=iturn3_end+2
iturn3_start=iturn3_start-1
iturn3_end=iturn3_end-1
+ call int_bounds(nres-3,itau_start,itau_end)
+ itau_start=itau_start+3
+ itau_end=itau_end+3
call int_bounds(nres-3,iphi1_start,iphi1_end)
iphi1_start=iphi1_start+3
iphi1_end=iphi1_end+3
idihconstr_end=ndih_constr
iphid_start=iphi_start
iphid_end=iphi_end-1
+ itau_start=4
+ itau_end=nres
ibond_start=2
ibond_end=nres-1
ibondp_start=nnt
include 'COMMON.INTERACT'
include 'COMMON.MD'
include 'COMMON.IOUNITS'
-
+ include 'COMMON.SCCOR'
c calculating dE/ddc1
- if (nres.lt.3) return
+ if (nres.lt.3) go to 18
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)
enddo
endif
enddo
+c----------------------------------------------------------------------
+C INTERTYP=1 SC...Ca...Ca...Ca
+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).and.(itype(1).ne.ntyp1)) then
+ do j=1,3
+cc Derviative was calculated for oposite vector of side chain therefore
+c there is "-" sign before gloc_sc
+ gxcart(j,1)=gxcart(j,1)-gloc_sc(1,0,icg)*
+ & 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).and.(itype(2).ne.ntyp1)) 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)
+ endif
+ enddo
+ endif
+ if ((nres.ge.3).and.(itype(3).ne.10).and.(itype(3).ne.ntyp1))
+ & then
+ do j=1,3
+ gcart(j,1)=gcart(j,1)+gloc_sc(2,1,icg)*dtauangle(j,2,1,4)
+ enddo
+ endif
+c As potetnial DO NOT depend on omicron anlge their derivative is
+c ommited
+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).and.(itype(2).ne.ntyp1)) 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).and.(itype(3).ne.ntyp1))
+ & 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(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(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).and.(itype(1).ne.ntyp1)) 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
+ gcart(j,2)=gcart(j,2)+gloc_sc(2,1,icg)*dtauangle(j,2,2,4)
+c 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).and.(itype(i).ne.ntyp1)) 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,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)
+ endif
+ if (itype(i+1).ne.10) then
+ gxcart(j,i)=gxcart(j,i)-gloc_sc(3,i-1,icg)
+ &*dtauangle(j,3,1,i+2)
+ gcart(j,i)=gcart(j,i)+gloc_sc(3,i-1,icg)
+ &*dtauangle(j,3,2,i+2)
+ endif
+ endif
+ if (itype(i-1).ne.10) then
+ gcart(j,i)=gcart(j,i)+gloc_sc(1,i-2,icg)*
+ & dtauangle(j,1,3,i+1)
+ endif
+ 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)*
+ & dtauangle(j,2,1,i+3)
+ endif
+ enddo
+ enddo
+ endif
+c Setting dE/ddnres-1
+ if(nres.ge.4) then
+ do j=1,3
+ if ((itype(nres-1).ne.10).and.(itype(nres-1).ne.ntyp1)) then
+ gxcart(j,nres-1)=gxcart(j,nres-1)+gloc_sc(2,nres-3,icg)
+ & *dtauangle(j,2,3,nres)
+c write (iout,*) "gxcart(nres-1)", gloc_sc(2,nres-3,icg),
+c & 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).and.(itype(nres).ne.ntyp1)) 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).and.(itype(nres-2).ne.ntyp1)) 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).and.(itype(nres).ne.ntyp1)) then
+ gcart(j,nres-1)=gcart(j,nres-1)+gloc_sc(2,nres-2,icg)*
+ & dtauangle(j,2,2,nres+1)
+c write (iout,*) "gcart(nres-1)", gloc_sc(2,nres-2,icg),
+c & 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).and.
+ & (itype(nres).ne.ntyp1))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
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),
if (itype(i-1).ne.21) dtheta(j,2,i)=-dcostheta(j,2,i)/sint
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).and.(itype(i-1).ne.ntyp1)) 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-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
+ 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)
+ & +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
+ enddo
+
c Derivatives of phi:
c If phi is 0 or 180 degrees, then the formulas
c have to be derived by power series expansion of the
enddo
endif
enddo
+Calculate derivative of Tauangle
+#ifdef PARINTDER
+ do i=itau_start,itau_end
+#else
+ do i=3,nres
+#endif
+ if ((itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10)) cycle
+c if ((itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10).or.
+c & (itype(i-1).eq.ntyp1).or.(itype(i).eq.ntyp1)) cycle
+cc dtauangle(j,intertyp,dervityp,residue number)
+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)
+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)
+ fac1=cost*fac0
+ 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
+ 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,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,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*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,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,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
+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
+ do i=itau_start,itau_end
+#else
+ do i=4,nres
+#endif
+ if ((itype(i-1).eq.ntyp1).or.(itype(i-1).eq.10).or.
+ & (itype(i-2).eq.ntyp1).or.(itype(i-3).eq.ntyp1)) 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))
+ cost1=dcos(theta(i-1))
+ cosg=dcos(tauangle(2,i))
+c do j=1,3
+c dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres)
+c enddo
+ scalp=scalar(dc_norm(1,i-3),dc_norm(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_norm(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)
+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)
+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-1+nres)
+c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
+ dtauangle(j,2,3,i)=cosg_inv*dsintau(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_norm(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)
+ dtauangle(j,2,2,i)=-1/sing*dcostau(j,2,2,i)
+ dcostau(j,2,3,i)=fac2*dcosomicron(j,1,2,i)+fac4*
+ & dcosomicron(j,1,2,i)-fac0*(dc_norm(j,i-3)-scalp*
+ & dc_norm(j,i-1+nres))/vbld(i-1+nres)
+ dtauangle(j,2,3,i)=-1/sing*dcostau(j,2,3,i)
+c write(iout,*) i,j,"else", dtauangle(j,2,3,i)
+ enddo
+ endif
+ enddo
+
+CCC third case SC...Ca...Ca...SC
+#ifdef PARINTDER
+
+ do i=itau_start,itau_end
+#else
+ do i=3,nres
+#endif
+c the conventional case
+ if ((itype(i-1).eq.ntyp1).or.(itype(i-1).eq.10).or.
+ &(itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10)) cycle
+ 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)
+c dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres)
+ enddo
+ scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(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_norm(1,i-1+nres),dc_norm(1,i-2),vp1)
+ call vecpr(dc_norm2(1,i-2+nres),dc_norm(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,2,i-1)
+ & -(fac0*vp1(j)-sing*dc_norm(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,1,i-1)+ctgt*domicron(j,1,1,i))
+ & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
+ 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_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,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,2,i-1)+fac3*
+ & dcosomicron(j,2,2,i-1)-fac0*(dc_norm(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,1,i-1)+fac2*
+ & dcosomicron(j,1,1,i)+fac3*dcosomicron(j,2,1,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*
+ & dcosomicron(j,1,2,i)-fac0*(dc_norm2(j,i-2+nres)-scalp*
+ & dc_norm(j,i-1+nres))/vbld(i-1+nres)
+ dtauangle(j,3,3,i)=-1/sing*dcostau(j,3,3,i)
+c write(iout,*) "else",i
+ enddo
+ endif
+ enddo
+
#ifdef CRYST_SC
c Derivatives of side-chain angles alpha and omega
#if defined(MPI) && defined(PARINTDER)
incr(j)=d_t(j,0)
enddo
do i=nnt,nct
- iti=itype(i)
+ iti=iabs(itype(i))
if (itype(i).eq.10) then
do j=1,3
v(j)=incr(j)
ind=ind+1
ii = ind+m
iti=itype(i)
- massvec(ii)=msc(iti)
- if (iti.ne.10 .and. iti.ne.21) then
+ massvec(ii)=msc(iabs(iti))
+ if (iti.ne.10 .and. iti.ne.ntyp1) then
ind1=ind1+1
ii1= ind1+m1
A(ii,ii1)=1.0d0
- Gmat(ii1,ii1)=ISC(iti)
+ Gmat(ii1,ii1)=ISC(iabs(iti))
endif
enddo
c Off-diagonal elements of the dX part of A
enddo
M_SC=0.0d0
do i=nnt,nct
- iti=itype(i)
- M_SC=M_SC+msc(iti)
+ iti=iabs(itype(i))
+ M_SC=M_SC+msc(iabs(iti))
inres=i+nres
do j=1,3
- cm(j)=cm(j)+msc(iti)*c(j,inres)
+ cm(j)=cm(j)+msc(iabs(iti))*c(j,inres)
enddo
enddo
do j=1,3
do j=1,3
pr(j)=c(j,inres)-cm(j)
enddo
- Im(1,1)=Im(1,1)+msc(iti)*(pr(2)*pr(2)+pr(3)*pr(3))
- Im(1,2)=Im(1,2)-msc(iti)*pr(1)*pr(2)
- Im(1,3)=Im(1,3)-msc(iti)*pr(1)*pr(3)
- Im(2,3)=Im(2,3)-msc(iti)*pr(2)*pr(3)
- Im(2,2)=Im(2,2)+msc(iti)*(pr(3)*pr(3)+pr(1)*pr(1))
- Im(3,3)=Im(3,3)+msc(iti)*(pr(1)*pr(1)+pr(2)*pr(2))
+ Im(1,1)=Im(1,1)+msc(iabs(iti))*(pr(2)*pr(2)+pr(3)*pr(3))
+ Im(1,2)=Im(1,2)-msc(iabs(iti))*pr(1)*pr(2)
+ Im(1,3)=Im(1,3)-msc(iabs(iti))*pr(1)*pr(3)
+ Im(2,3)=Im(2,3)-msc(iabs(iti))*pr(2)*pr(3)
+ Im(2,2)=Im(2,2)+msc(iabs(iti))*(pr(3)*pr(3)+pr(1)*pr(1))
+ Im(3,3)=Im(3,3)+msc(iabs(iti))*(pr(1)*pr(1)+pr(2)*pr(2))
enddo
do i=nnt,nct-1
c write (iout,*) "i",i," iti",iti," pr",(pr(j),j=1,3),
c & " v",(v(j),j=1,3)," vp",(vp(j),j=1,3)
do j=1,3
- L(j)=L(j)+msc(iti)*vp(j)
+ L(j)=L(j)+msc(iabs(iti))*vp(j)
enddo
c write (iout,*) "L",(l(j),j=1,3)
if (itype(i).ne.10 .and. itype(i).ne.21) then
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=119,err=119) nsccortyp
+#ifdef SCCORPDB
+ read (isccor,*,end=119,err=119) (isccortyp(i),i=1,ntyp)
+ do i=-ntyp,-1
+ isccortyp(i)=-isccortyp(-i)
+ enddo
+ iscprol=isccortyp(20)
+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=119,err=119)
+ &nterm_sccor(i,j),nlor_sccor(i,j)
+ v0ijsccor=0.0d0
+ si=-1.0d0
+ nterm_sccor(-i,j)=nterm_sccor(i,j)
+ nterm_sccor(-i,-j)=nterm_sccor(i,j)
+ nterm_sccor(i,-j)=nterm_sccor(i,j)
+ do k=1,nterm_sccor(i,j)
+ read (isccor,*,end=119,err=119) kk,v1sccor(k,l,i,j)
+ & ,v2sccor(k,l,i,j)
+ if (j.eq.iscprol) then
+ if (i.eq.isccortyp(10)) then
+ v1sccor(k,l,i,-j)=v1sccor(k,l,i,j)
+ v2sccor(k,l,i,-j)=-v2sccor(k,l,i,j)
+ else
+ v1sccor(k,l,i,-j)=v1sccor(k,l,i,j)*0.5d0
+ & +v2sccor(k,l,i,j)*dsqrt(0.75d0)
+ v2sccor(k,l,i,-j)=-v2sccor(k,l,i,j)*0.5d0
+ & +v1sccor(k,l,i,j)*dsqrt(0.75d0)
+ v1sccor(k,l,-i,-j)=v1sccor(k,l,i,j)
+ v2sccor(k,l,-i,-j)=-v2sccor(k,l,i,j)
+ v1sccor(k,l,-i,j)=v1sccor(k,l,i,-j)
+ v2sccor(k,l,-i,j)=-v2sccor(k,l,i,-j)
+ endif
+ else
+ if (i.eq.isccortyp(10)) then
+ v1sccor(k,l,i,-j)=v1sccor(k,l,i,j)
+ v2sccor(k,l,i,-j)=-v2sccor(k,l,i,j)
+ else
+ if (j.eq.isccortyp(10)) then
+ v1sccor(k,l,-i,j)=v1sccor(k,l,i,j)
+ v2sccor(k,l,-i,j)=-v2sccor(k,l,i,j)
+ else
+ v1sccor(k,l,i,-j)=-v1sccor(k,l,i,j)
+ v2sccor(k,l,i,-j)=-v2sccor(k,l,i,j)
+ v1sccor(k,l,-i,-j)=v1sccor(k,l,i,j)
+ v2sccor(k,l,-i,-j)=-v2sccor(k,l,i,j)
+ v1sccor(k,l,-i,j)=v1sccor(k,l,i,-j)
+ v2sccor(k,l,-i,j)=-v2sccor(k,l,i,-j)
+ endif
+ endif
+ endif
+ v0ijsccor=v0ijsccor+si*v1sccor(k,l,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)
+#else
+ 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,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)
+ 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)
+
+#endif
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,'(2(1pe15.5))') v1sccor(k,i,j),v2sccor(k,i,j)
+ write (iout,*) 'Fourier constants'
+ do k=1,nterm_sccor(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)
+ write (iout,'(3(1pe15.5))')
+ & vlor1sccor(k,i,j),vlor2sccor(k,i,j),vlor3sccor(k,i,j)
enddo
enddo
enddo
endif
+
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.
ind=ind+3
enddo
do i=nnt,nct
+<<<<<<< HEAD
if (itype(i).ne.10 .and. itype(i).ne.21) then
+=======
+ if (itype(i).ne.10) then
+>>>>>>> 2acc991... Wprowadznie potencjalow SC-COR do multichain oraz ich pseudosymetrii dla
do j=1,3
d_t_work(ind+j)=d_t(j,i+nres)
enddo
ind=ind+3
enddo
do i=nnt,nct
+<<<<<<< HEAD
if (itype(i).ne.10 .and. itype(i).ne.21) then
+=======
+ if (itype(i).ne.10) then
+>>>>>>> 2acc991... Wprowadznie potencjalow SC-COR do multichain oraz ich pseudosymetrii dla
do j=1,3
friction(j,i+nres)=fric_work(ind+j)
enddo
stochforc(j,0)=ff(j)+force(j,nnt+nres)
enddo
do i=nnt,nct
+<<<<<<< HEAD
if (itype(i).ne.10 .and. itype(i).ne.21) then
+=======
+ if (itype(i).ne.10) then
+>>>>>>> 2acc991... Wprowadznie potencjalow SC-COR do multichain oraz ich pseudosymetrii dla
do j=1,3
stochforc(j,i+nres)=force(j,i+nres)
enddo
ind=ind+3
enddo
do i=nnt,nct
+<<<<<<< HEAD
if (itype(i).ne.10 .and. itype(i).ne.21) then
+=======
+ if (itype(i).ne.10) then
+>>>>>>> 2acc991... Wprowadznie potencjalow SC-COR do multichain oraz ich pseudosymetrii dla
do j=1,3
stochforcvec(ind+j)=stochforc(j,i+nres)
enddo
c Load the friction coefficients corresponding to side chains
m=nct-nnt
ind=0
- gamsc(ntyp1)=1.0d0
do i=nnt,nct
ind=ind+1
ii = ind+m
time00=tcpu()
#endif
call MPI_Bcast(10,1,MPI_INTEGER,king,FG_COMM,IERROR)
+#ifdef MPI
time_Bcast=time_Bcast+MPI_Wtime()-time00
+#else
+ time_Bcast=time_Bcast+tcpu()-time00
+#endif
c print *,"Processor",myrank,
c & " BROADCAST iorder in SETUP_FRICMAT"
endif
c licznik=licznik+1
c write (iout,*) "setup_fricmat licznik",licznik
+#ifdef MPI
time00=MPI_Wtime()
+#else
+ time00=tcpu()
+#endif
c Scatter the friction matrix
call MPI_Scatterv(fricmat(1,1),nginv_counts(0),
& nginv_start(0),MPI_DOUBLE_PRECISION,fcopy(1,1),
& myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
- time_scatter=time_scatter+MPI_Wtime()-time00
#ifdef TIMING
#ifdef MPI
time_scatter_fmat=time_scatter_fmat+MPI_Wtime()-time00
#else
+ time_scatter=time_scatter+tcpu()-time00
time_scatter_fmat=time_scatter_fmat+tcpu()-time00
#endif
#endif
include 'COMMON.NAMES'
double precision radius(maxres2),gamvec(maxres2)
parameter (twosix=1.122462048309372981d0)
- logical lprn /.true./
+ logical lprn /.false./
c
c determine new friction coefficients every few SD steps
c
& dcostau,dsintau,dtauangle,dcosomicron,
& domicron,v0sccor
integer nterm_sccor,isccortyp,nsccortyp,nlor_sccor
+<<<<<<< HEAD
common/sccor/v1sccor(maxterm_sccor,3,20,20),
& v2sccor(maxterm_sccor,3,20,20),
+=======
+ common/sccor/v1sccor(maxterm_sccor,3,-ntyp:ntyp,-ntyp:ntyp),
+ & v2sccor(maxterm_sccor,3,-ntyp:ntyp,-ntyp:ntyp),
+ & v0sccor(maxterm_sccor,-ntyp:ntyp),
+ & nterm_sccor(-ntyp:ntyp,-ntyp:ntyp),isccortyp(-ntyp:ntyp),
+ & nsccortyp,
+ & nlor_sccor(-ntyp:ntyp,-ntyp:ntyp),
+>>>>>>> 2acc991... Wprowadznie potencjalow SC-COR do multichain oraz ich pseudosymetrii dla
& vlor1sccor(maxterm_sccor,20,20),
& vlor2sccor(maxterm_sccor,20,20),
& vlor3sccor(maxterm_sccor,20,20),gloc_sc(3,0:maxres2,10),
projects / unres.git / commitdiff