projects
/
unres.git
/ commitdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
| commitdiff |
tree
raw
|
patch
|
inline
| side by side (parent:
02826f7
)
dzialajacy gradient
author
Adam Sieradzan
<adasko@piasek4.chem.univ.gda.pl>
Mon, 18 Nov 2013 19:02:48 +0000
(20:02 +0100)
committer
Adam Sieradzan
<adasko@piasek4.chem.univ.gda.pl>
Mon, 18 Nov 2013 19:02:48 +0000
(20:02 +0100)
source/unres/src_MD-M/CMakeLists.txt
patch
|
blob
|
history
source/unres/src_MD-M/COMMON.CONTACTS
patch
|
blob
|
history
source/unres/src_MD-M/COMMON.TORSION
patch
|
blob
|
history
source/unres/src_MD-M/energy_p_new_barrier.F
patch
|
blob
|
history
source/unres/src_MD-M/parmread.F
patch
|
blob
|
history
diff --git
a/source/unres/src_MD-M/CMakeLists.txt
b/source/unres/src_MD-M/CMakeLists.txt
index
cfdc654
..
3884bc1
100644
(file)
--- a/
source/unres/src_MD-M/CMakeLists.txt
+++ b/
source/unres/src_MD-M/CMakeLists.txt
@@
-203,7
+203,7
@@
set_property(SOURCE ${UNRES_MDM_SRC3} PROPERTY COMPILE_FLAGS ${FFLAGS3} )
#=========================================
if(UNRES_MD_FF STREQUAL "GAB" )
# set preprocesor flags
#=========================================
if(UNRES_MD_FF STREQUAL "GAB" )
# set preprocesor flags
- set(CPPFLAGS "PROCOR -DUNRES -DISNAN -DSPLITELE -DLANG0 -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORPDB" )
+ set(CPPFLAGS "PROCOR -DUNRES -DISNAN -DSPLITELE -DLANG0 -DCRYST_BOND -DCRYST_THETA -DCRYST_SC -DSCCORPDB -DNEWCORR" )
#=========================================
# Settings for E0LL2Y force field
#=========================================
# Settings for E0LL2Y force field
diff --git
a/source/unres/src_MD-M/COMMON.CONTACTS
b/source/unres/src_MD-M/COMMON.CONTACTS
index
4b81714
..
45c578b
100644
(file)
--- a/
source/unres/src_MD-M/COMMON.CONTACTS
+++ b/
source/unres/src_MD-M/COMMON.CONTACTS
@@
-28,7
+28,8
@@
C 10/30/99 Added other pre-computed vectors and matrices needed
C to calculate three - six-order el-loc correlation terms
double precision Ug,Ugder,Ug2,Ug2der,obrot,obrot2,obrot_der,
& obrot2_der,Ub2,Ub2der,mu,muder,EUg,EUgder,CUg,CUgder,gmu,gUb2
C to calculate three - six-order el-loc correlation terms
double precision Ug,Ugder,Ug2,Ug2der,obrot,obrot2,obrot_der,
& obrot2_der,Ub2,Ub2der,mu,muder,EUg,EUgder,CUg,CUgder,gmu,gUb2
- & DUg,DUgder,DtUg2,DtUg2der,Ctobr,Ctobrder,Dtobr2,Dtobr2der
+ & DUg,DUgder,DtUg2,DtUg2der,Ctobr,Ctobrder,Dtobr2,Dtobr2der,
+ & gtEug
common /rotat/ Ug(2,2,maxres),Ugder(2,2,maxres),Ug2(2,2,maxres),
& Ug2der(2,2,maxres),obrot(2,maxres),obrot2(2,maxres),
& obrot_der(2,maxres),obrot2_der(2,maxres)
common /rotat/ Ug(2,2,maxres),Ugder(2,2,maxres),Ug2(2,2,maxres),
& Ug2der(2,2,maxres),obrot(2,maxres),obrot2(2,maxres),
& obrot_der(2,maxres),obrot2_der(2,maxres)
@@
-40,7
+41,7
@@
C amino-acid residue.
& Dtobr2(2,maxres),Dtobr2der(2,maxres),
& EUg(2,2,maxres),EUgder(2,2,maxres),CUg(2,2,maxres),
& CUgder(2,2,maxres),DUg(2,2,maxres),Dugder(2,2,maxres),
& Dtobr2(2,maxres),Dtobr2der(2,maxres),
& EUg(2,2,maxres),EUgder(2,2,maxres),CUg(2,2,maxres),
& CUgder(2,2,maxres),DUg(2,2,maxres),Dugder(2,2,maxres),
- & DtUg2(2,2,maxres),DtUg2der(2,2,maxres)
+ & DtUg2(2,2,maxres),DtUg2der(2,2,maxres),gtEUg(2,2,maxres)
C This common block contains vectors and matrices dependent on two
C consecutive amino-acid residues.
double precision Ug2Db1t,Ug2Db1tder,CUgb2,CUgb2der,EUgC,
C This common block contains vectors and matrices dependent on two
C consecutive amino-acid residues.
double precision Ug2Db1t,Ug2Db1tder,CUgb2,CUgb2der,EUgC,
diff --git
a/source/unres/src_MD-M/COMMON.TORSION
b/source/unres/src_MD-M/COMMON.TORSION
index
95472be
..
3f1b981
100644
(file)
--- a/
source/unres/src_MD-M/COMMON.TORSION
+++ b/
source/unres/src_MD-M/COMMON.TORSION
@@
-25,16
+25,17
@@
C 6/23/01 - constants for double torsionals
C 9/18/99 - added Fourier coeffficients of the expansion of local energy
C surface
double precision b1,b2,cc,dd,ee,ctilde,dtilde,b2tilde,b1tilde,
C 9/18/99 - added Fourier coeffficients of the expansion of local energy
C surface
double precision b1,b2,cc,dd,ee,ctilde,dtilde,b2tilde,b1tilde,
- &bnew1,bnew2,eenew,ggb1,ggb2
+ &bnew1,bnew2,eenew,gtb1,gtb2,eeold,gtee
integer nloctyp
common/fourier/ b1(2,maxres),b2(2,maxres),
& bnew1(3,2,-maxtor:maxtor),bnew2(3,2,-maxtor:maxtor),
& cc(2,2,-maxtor:maxtor),
integer nloctyp
common/fourier/ b1(2,maxres),b2(2,maxres),
& bnew1(3,2,-maxtor:maxtor),bnew2(3,2,-maxtor:maxtor),
& cc(2,2,-maxtor:maxtor),
- & dd(2,2,-maxtor:maxtor),ee(2,2,-maxtor:maxtor),
+ & dd(2,2,-maxtor:maxtor),eeold(2,2,-maxtor:maxtor),
& eenew(2,-maxtor:maxtor),
& eenew(2,-maxtor:maxtor),
+ & ee(2,2,maxres),
& ctilde(2,2,-maxtor:maxtor),
& dtilde(2,2,-maxtor:maxtor),b1tilde(2,maxres),
& b2tilde(2,maxres),
& ctilde(2,2,-maxtor:maxtor),
& dtilde(2,2,-maxtor:maxtor),b1tilde(2,maxres),
& b2tilde(2,maxres),
- & gtb1(2,maxres),gtb2(2,maxres),
+ & gtb1(2,maxres),gtb2(2,maxres),gtEE(2,2,maxres),
& nloctyp
& 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
2682840
..
27c5f1d
100644
(file)
--- a/
source/unres/src_MD-M/energy_p_new_barrier.F
+++ b/
source/unres/src_MD-M/energy_p_new_barrier.F
@@
-2297,7
+2297,14
@@
c endif
gtb1(2,i-2)=0.0
b2(2,i-2)=bnew2(1,2,iti)
gtb2(2,i-2)=0.0
gtb1(2,i-2)=0.0
b2(2,i-2)=bnew2(1,2,iti)
gtb2(2,i-2)=0.0
-c EE(1,1,iti)=0.0d0
+ EE(1,1,i-2)=eenew(1,iti)*dcos(theta(i-1))
+ EE(1,2,i-2)=eeold(1,2,iti)
+ EE(2,1,i-2)=eeold(2,1,iti)
+ EE(2,2,i-2)=eeold(2,2,iti)
+ gtEE(1,1,i-2)=-eenew(1,iti)*dsin(theta(i-1))
+ gtEE(1,2,i-2)=0.0d0
+ gtEE(2,2,i-2)=0.0d0
+ gtEE(2,1,i-2)=0.0d0
c EE(2,2,iti)=0.0d0
c EE(1,2,iti)=0.5d0*eenew(1,iti)
c EE(2,1,iti)=0.5d0*eenew(1,iti)
c EE(2,2,iti)=0.0d0
c EE(1,2,iti)=0.5d0*eenew(1,iti)
c EE(2,1,iti)=0.5d0*eenew(1,iti)
@@
-2384,7
+2391,6
@@
c write (iout,*) 'theta=', theta(i-1)
Ug2der(2,2,i-2)=0.0d0
endif
c if (i.gt. iatel_s+2 .and. i.lt.iatel_e+5) then
Ug2der(2,2,i-2)=0.0d0
endif
c if (i.gt. iatel_s+2 .and. i.lt.iatel_e+5) then
-#ifndef NEWCORR
if (i.gt. nnt+2 .and. i.lt.nct+2) then
iti = itortyp(itype(i-2))
else
if (i.gt. nnt+2 .and. i.lt.nct+2) then
iti = itortyp(itype(i-2))
else
@@
-2396,7
+2402,6
@@
c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then
else
iti1=ntortyp+1
endif
else
iti1=ntortyp+1
endif
-#endif
cd write (iout,*) '*******i',i,' iti1',iti
cd write (iout,*) 'b1',b1(:,iti)
cd write (iout,*) 'b2',b2(:,iti)
cd write (iout,*) '*******i',i,' iti1',iti
cd write (iout,*) 'b1',b1(:,iti)
cd write (iout,*) 'b2',b2(:,iti)
@@
-2408,7
+2413,13
@@
c if (i .gt. iatel_s+2) then
call matvec2(Ug(1,1,i-2),gtb2(1,i-2),gUb2(1,i-2))
c write (iout,*) Ug(1,1,i-2),gtb2(1,i-2),gUb2(1,i-2),"chuj"
#endif
call matvec2(Ug(1,1,i-2),gtb2(1,i-2),gUb2(1,i-2))
c write (iout,*) Ug(1,1,i-2),gtb2(1,i-2),gUb2(1,i-2),"chuj"
#endif
- call matmat2(EE(1,1,iti),Ug(1,1,i-2),EUg(1,1,i-2))
+c write(iout,*) "co jest kurwa", iti, EE(1,1,iti),EE(2,1,iti),
+c & EE(1,2,iti),EE(2,2,iti)
+ call matmat2(EE(1,1,i-2),Ug(1,1,i-2),EUg(1,1,i-2))
+ call matmat2(gtEE(1,1,i-2),Ug(1,1,i-2),gtEUg(1,1,i-2))
+c write(iout,*) "Macierz EUG",
+c & eug(1,1,i-2),eug(1,2,i-2),eug(2,1,i-2),
+c & eug(2,2,i-2)
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0)
& then
call matmat2(CC(1,1,iti),Ug(1,1,i-2),CUg(1,1,i-2))
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0)
& then
call matmat2(CC(1,1,iti),Ug(1,1,i-2),CUg(1,1,i-2))
@@
-2431,7
+2442,7
@@
c write (iout,*) Ug(1,1,i-2),gtb2(1,i-2),gUb2(1,i-2),"chuj"
enddo
endif
call matvec2(Ugder(1,1,i-2),b2(1,i-2),Ub2der(1,i-2))
enddo
endif
call matvec2(Ugder(1,1,i-2),b2(1,i-2),Ub2der(1,i-2))
- call matmat2(EE(1,1,iti),Ugder(1,1,i-2),EUgder(1,1,i-2))
+ call matmat2(EE(1,1,i-2),Ugder(1,1,i-2),EUgder(1,1,i-2))
do k=1,2
muder(k,i-2)=Ub2der(k,i-2)
enddo
do k=1,2
muder(k,i-2)=Ub2der(k,i-2)
enddo
@@
-2867,8
+2878,7
@@
C
ymedi=c(2,i)+0.5d0*dyi
zmedi=c(3,i)+0.5d0*dzi
num_conti=0
ymedi=c(2,i)+0.5d0*dyi
zmedi=c(3,i)+0.5d0*dzi
num_conti=0
-c TU ZLE
-c call eelecij(i,i+2,ees,evdw1,eel_loc)
+ call eelecij(i,i+2,ees,evdw1,eel_loc)
if (wturn3.gt.0.0d0) call eturn3(i,eello_turn3)
num_cont_hb(i)=num_conti
enddo
if (wturn3.gt.0.0d0) call eturn3(i,eello_turn3)
num_cont_hb(i)=num_conti
enddo
@@
-2886,8
+2896,8
@@
c call eelecij(i,i+2,ees,evdw1,eel_loc)
ymedi=c(2,i)+0.5d0*dyi
zmedi=c(3,i)+0.5d0*dzi
num_conti=num_cont_hb(i)
ymedi=c(2,i)+0.5d0*dyi
zmedi=c(3,i)+0.5d0*dzi
num_conti=num_cont_hb(i)
-c TU ZLE
-c call eelecij(i,i+3,ees,evdw1,eel_loc)
+c write(iout,*) "JESTEM W PETLI"
+ call eelecij(i,i+3,ees,evdw1,eel_loc)
if (wturn4.gt.0.0d0 .and. itype(i+2).ne.ntyp1)
& call eturn4(i,eello_turn4)
num_cont_hb(i)=num_conti
if (wturn4.gt.0.0d0 .and. itype(i+2).ne.ntyp1)
& call eturn4(i,eello_turn4)
num_cont_hb(i)=num_conti
@@
-2895,8
+2905,8
@@
c call eelecij(i,i+3,ees,evdw1,eel_loc)
c
c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
c
c
c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
c
-c do i=iatel_s,iatel_e
- do i=7,7
+ do i=iatel_s,iatel_e
+c do i=7,7
if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle
dxi=dc(1,i)
dyi=dc(2,i)
if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle
dxi=dc(1,i)
dyi=dc(2,i)
@@
-2909,8
+2919,8
@@
c do i=iatel_s,iatel_e
zmedi=c(3,i)+0.5d0*dzi
c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
num_conti=num_cont_hb(i)
zmedi=c(3,i)+0.5d0*dzi
c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
num_conti=num_cont_hb(i)
-c do j=ielstart(i),ielend(i)
- do j=13,13
+ do j=ielstart(i),ielend(i)
+c do j=13,13
c write (iout,*) 'tu wchodze',i,j,itype(i),itype(j)
if (itype(j).eq.ntyp1.or. itype(j+1).eq.ntyp1) cycle
call eelecij(i,j,ees,evdw1,eel_loc)
c write (iout,*) 'tu wchodze',i,j,itype(i),itype(j)
if (itype(j).eq.ntyp1.or. itype(j+1).eq.ntyp1) cycle
call eelecij(i,j,ees,evdw1,eel_loc)
@@
-3652,7
+3662,9
@@
C Third- and fourth-order contributions from turns
dimension ggg(3)
double precision auxmat(2,2),auxmat1(2,2),auxmat2(2,2),pizda(2,2),
& e1t(2,2),e2t(2,2),e3t(2,2),e1tder(2,2),e2tder(2,2),e3tder(2,2),
dimension ggg(3)
double precision auxmat(2,2),auxmat1(2,2),auxmat2(2,2),pizda(2,2),
& e1t(2,2),e2t(2,2),e3t(2,2),e1tder(2,2),e2tder(2,2),e3tder(2,2),
- & e1a(2,2),ae3(2,2),ae3e2(2,2),auxvec(2),auxvec1(2)
+ & e1a(2,2),ae3(2,2),ae3e2(2,2),auxvec(2),auxvec1(2),gpizda1(2,2),
+ & gpizda2(2,2),auxgmat1(2,2),auxgmatt1(2,2),
+ & auxgmat2(2,2),auxgmatt2(2,2)
double precision agg(3,4),aggi(3,4),aggi1(3,4),
& aggj(3,4),aggj1(3,4),a_temp(2,2),auxmat3(2,2)
common /locel/ a_temp,agg,aggi,aggi1,aggj,aggj1,a22,a23,a32,a33,
double precision agg(3,4),aggi(3,4),aggi1(3,4),
& aggj(3,4),aggj1(3,4),a_temp(2,2),auxmat3(2,2)
common /locel/ a_temp,agg,aggi,aggi1,aggj,aggj1,a22,a23,a32,a33,
@@
-3676,9
+3688,24
@@
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
cd call checkint_turn3(i,a_temp,eello_turn3_num)
call matmat2(EUg(1,1,i+1),EUg(1,1,i+2),auxmat(1,1))
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
cd call checkint_turn3(i,a_temp,eello_turn3_num)
call matmat2(EUg(1,1,i+1),EUg(1,1,i+2),auxmat(1,1))
+c auxalary matices for theta gradient
+c auxalary matrix for i+1 and constant i+2
+ call matmat2(gtEUg(1,1,i+1),EUg(1,1,i+2),auxgmat1(1,1))
+c auxalary matrix for i+2 and constant i+1
+ call matmat2(EUg(1,1,i+1),gtEUg(1,1,i+2),auxgmat2(1,1))
call transpose2(auxmat(1,1),auxmat1(1,1))
call transpose2(auxmat(1,1),auxmat1(1,1))
+ call transpose2(auxgmat1(1,1),auxgmatt1(1,1))
+ call transpose2(auxgmat2(1,1),auxgmatt2(1,1))
call matmat2(a_temp(1,1),auxmat1(1,1),pizda(1,1))
call matmat2(a_temp(1,1),auxmat1(1,1),pizda(1,1))
+ call matmat2(a_temp(1,1),auxgmatt1(1,1),gpizda1(1,1))
+ call matmat2(a_temp(1,1),auxgmatt2(1,1),gpizda2(1,1))
eello_turn3=eello_turn3+0.5d0*(pizda(1,1)+pizda(2,2))
eello_turn3=eello_turn3+0.5d0*(pizda(1,1)+pizda(2,2))
+C Derivatives in theta
+ gloc(nphi+i,icg)=gloc(nphi+i,icg)
+ & +0.5d0*(gpizda1(1,1)+gpizda1(2,2))*wturn3
+ gloc(nphi+i+1,icg)=gloc(nphi+i+1,icg)
+ & +0.5d0*(gpizda2(1,1)+gpizda2(2,2))*wturn3
+
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
& 'eturn3',i,j,0.5d0*(pizda(1,1)+pizda(2,2))
cd write (2,*) 'i,',i,' j',j,'eello_turn3',
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
& 'eturn3',i,j,0.5d0*(pizda(1,1)+pizda(2,2))
cd write (2,*) 'i,',i,' j',j,'eello_turn3',
@@
-3752,7
+3779,11
@@
C Third- and fourth-order contributions from turns
dimension ggg(3)
double precision auxmat(2,2),auxmat1(2,2),auxmat2(2,2),pizda(2,2),
& e1t(2,2),e2t(2,2),e3t(2,2),e1tder(2,2),e2tder(2,2),e3tder(2,2),
dimension ggg(3)
double precision auxmat(2,2),auxmat1(2,2),auxmat2(2,2),pizda(2,2),
& e1t(2,2),e2t(2,2),e3t(2,2),e1tder(2,2),e2tder(2,2),e3tder(2,2),
- & e1a(2,2),ae3(2,2),ae3e2(2,2),auxvec(2),auxvec1(2),auxgvec(2)
+ & e1a(2,2),ae3(2,2),ae3e2(2,2),auxvec(2),auxvec1(2),auxgvec(2),
+ & auxgEvec1(2),auxgEvec2(2),auxgEvec3(2),
+ & gte1t(2,2),gte2t(2,2),gte3t(2,2),
+ & gte1a(2,2),gtae3(2,2),gtae3e2(2,2), ae3gte2(2,2),
+ & gtEpizda1(2,2),gtEpizda2(2,2),gtEpizda3(2,2)
double precision agg(3,4),aggi(3,4),aggi1(3,4),
& aggj(3,4),aggj1(3,4),a_temp(2,2),auxmat3(2,2)
common /locel/ a_temp,agg,aggi,aggi1,aggj,aggj1,a22,a23,a32,a33,
double precision agg(3,4),aggi(3,4),aggi1(3,4),
& aggj(3,4),aggj1(3,4),a_temp(2,2),auxmat3(2,2)
common /locel/ a_temp,agg,aggi,aggi1,aggj,aggj1,a22,a23,a32,a33,
@@
-3772,6
+3803,7
@@
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
cd call checkint_turn4(i,a_temp,eello_turn4_num)
c write (iout,*) "eturn4 i",i," j",j," j1",j1," j2",j2
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
cd call checkint_turn4(i,a_temp,eello_turn4_num)
c write (iout,*) "eturn4 i",i," j",j," j1",j1," j2",j2
+c write(iout,*)"WCHODZE W PROGRAM"
a_temp(1,1)=a22
a_temp(1,2)=a23
a_temp(2,1)=a32
a_temp(1,1)=a22
a_temp(1,2)=a23
a_temp(2,1)=a32
@@
-3783,37
+3815,83
@@
c write(iout,*) "iti1",iti1," iti2",iti2," iti3",iti3
call transpose2(EUg(1,1,i+1),e1t(1,1))
call transpose2(Eug(1,1,i+2),e2t(1,1))
call transpose2(Eug(1,1,i+3),e3t(1,1))
call transpose2(EUg(1,1,i+1),e1t(1,1))
call transpose2(Eug(1,1,i+2),e2t(1,1))
call transpose2(Eug(1,1,i+3),e3t(1,1))
+C Ematrix derivative in theta
+ call transpose2(gtEUg(1,1,i+1),gte1t(1,1))
+ call transpose2(gtEug(1,1,i+2),gte2t(1,1))
+ call transpose2(gtEug(1,1,i+3),gte3t(1,1))
call matmat2(e1t(1,1),a_temp(1,1),e1a(1,1))
call matmat2(e1t(1,1),a_temp(1,1),e1a(1,1))
+c eta1 in derivative theta
+ call matmat2(gte1t(1,1),a_temp(1,1),gte1a(1,1))
call matvec2(e1a(1,1),Ub2(1,i+3),auxvec(1))
call matvec2(e1a(1,1),Ub2(1,i+3),auxvec(1))
+c auxgvec is derivative of Ub2 so i+3 theta
call matvec2(e1a(1,1),gUb2(1,i+3),auxgvec(1))
call matvec2(e1a(1,1),gUb2(1,i+3),auxgvec(1))
+c auxalary matrix of E i+1
+ call matvec2(gte1a(1,1),Ub2(1,i+3),auxgEvec1(1))
c s1=0.0
c gs1=0.0
s1=scalar2(b1(1,i+2),auxvec(1))
c s1=0.0
c gs1=0.0
s1=scalar2(b1(1,i+2),auxvec(1))
-c gs1=scalar2(gtb1(1,i+2),auxgvec(1))
+c derivative of theta i+2 with constant i+3
+ gs23=scalar2(gtb1(1,i+2),auxvec(1))
+c derivative of theta i+2 with constant i+2
+ gs32=scalar2(b1(1,i+2),auxgvec(1))
+c derivative of E matix in theta of i+1
+ gsE13=scalar2(b1(1,i+2),auxgEvec1(1))
+
call matmat2(a_temp(1,1),e3t(1,1),ae3(1,1))
call matmat2(a_temp(1,1),e3t(1,1),ae3(1,1))
+c ea31 in derivative theta
+ call matmat2(a_temp(1,1),gte3t(1,1),gtae3(1,1))
call matvec2(ae3(1,1),Ub2(1,i+2),auxvec(1))
call matvec2(ae3(1,1),Ub2(1,i+2),auxvec(1))
+c auxilary matrix auxgvec of Ub2 with constant E matirx
call matvec2(ae3(1,1),gUb2(1,i+2),auxgvec(1))
call matvec2(ae3(1,1),gUb2(1,i+2),auxgvec(1))
+c auxilary matrix auxgEvec1 of E matix with Ub2 constant
+ call matvec2(gtae3(1,1),Ub2(1,i+2),auxgEvec3(1))
+
c s2=0.0
c gs2=0.0
s2=scalar2(b1(1,i+1),auxvec(1))
c s2=0.0
c gs2=0.0
s2=scalar2(b1(1,i+1),auxvec(1))
-c gs2=scalar2(gtb1(1,i+1),auxgvec(1))
-c write(iout,*) gs1,gs2,'i=',i,auxgvec(1),gUb2(1,i+2),ggb1(1,i+2),
-c & ggb1(1,i+1)
+c derivative of theta i+1 with constant i+3
+ gs13=scalar2(gtb1(1,i+1),auxvec(1))
+c derivative of theta i+2 with constant i+1
+ gs21=scalar2(b1(1,i+1),auxgvec(1))
+c derivative of theta i+3 with constant i+1
+ gsE31=scalar2(b1(1,i+1),auxgEvec3(1))
+c write(iout,*) gs1,gs2,'i=',i,auxgvec(1),gUb2(1,i+2),gtb1(1,i+2),
+c & gtb1(1,i+1)
call matmat2(ae3(1,1),e2t(1,1),ae3e2(1,1))
call matmat2(ae3(1,1),e2t(1,1),ae3e2(1,1))
+c two derivatives over diffetent matrices
+c gtae3e2 is derivative over i+3
+ call matmat2(gtae3(1,1),e2t(1,1),gtae3e2(1,1))
+c ae3gte2 is derivative over i+2
+ call matmat2(ae3(1,1),gte2t(1,1),ae3gte2(1,1))
call matmat2(ae3e2(1,1),e1t(1,1),pizda(1,1))
call matmat2(ae3e2(1,1),e1t(1,1),pizda(1,1))
+c three possible derivative over theta E matices
+c i+1
+ call matmat2(ae3e2(1,1),gte1t(1,1),gtEpizda1(1,1))
+c i+2
+ call matmat2(ae3gte2(1,1),e1t(1,1),gtEpizda2(1,1))
+c i+3
+ call matmat2(gtae3e2(1,1),e1t(1,1),gtEpizda3(1,1))
s3=0.5d0*(pizda(1,1)+pizda(2,2))
s3=0.5d0*(pizda(1,1)+pizda(2,2))
+
+ gsEE1=0.5d0*(gtEpizda1(1,1)+gtEpizda1(2,2))
+ gsEE2=0.5d0*(gtEpizda2(1,1)+gtEpizda2(2,2))
+ gsEE3=0.5d0*(gtEpizda3(1,1)+gtEpizda3(2,2))
+
eello_turn4=eello_turn4-(s1+s2+s3)
#ifdef NEWCORR
eello_turn4=eello_turn4-(s1+s2+s3)
#ifdef NEWCORR
-c geel_loc_ij=-(gs1+gs2)
-c gloc(nphi+i,icg)=gloc(nphi+i,icg)-
-c & gs1
+ gloc(nphi+i,icg)=gloc(nphi+i,icg)
+ & -(gs13+gsE13+gsEE1)*wturn4
+ gloc(nphi+i+1,icg)= gloc(nphi+i+1,icg)
+ & -(gs23+gs21+gsEE2)*wturn4
+ gloc(nphi+i+2,icg)= gloc(nphi+i+2,icg)
+ & -(gs32+gsE31+gsEE3)*wturn4
c gloc(nphi+i+1,icg)=gloc(nphi+i+1,icg)-
c & gs2
#endif
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
& 'eturn4',i,j,-(s1+s2+s3)
c gloc(nphi+i+1,icg)=gloc(nphi+i+1,icg)-
c & gs2
#endif
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
& 'eturn4',i,j,-(s1+s2+s3)
-cd write (2,*) 'i,',i,' j',j,'eello_turn4',-(s1+s2+s3),
-cd & ' eello_turn4_num',8*eello_turn4_num
+c write (iout,*) 'i,',i,' j',j,'eello_turn4',-(s1+s2+s3),
+c & ' eello_turn4_num',8*eello_turn4_num
C Derivatives in gamma(i)
call transpose2(EUgder(1,1,i+1),e1tder(1,1))
call matmat2(e1tder(1,1),a_temp(1,1),auxmat(1,1))
C Derivatives in gamma(i)
call transpose2(EUgder(1,1,i+1),e1tder(1,1))
call matmat2(e1tder(1,1),a_temp(1,1),auxmat(1,1))
@@
-4556,7
+4634,7
@@
C Derivatives of the "mean" values in gamma1 and gamma2.
& 'ebend',i,ethetai
if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*E_tc*dthetg1
if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*E_tc*dthetg2
& 'ebend',i,ethetai
if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*E_tc*dthetg1
if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*E_tc*dthetg2
- gloc(nphi+i-2,icg)=wang*(E_theta+E_tc*dthett)
+ gloc(nphi+i-2,icg)=wang*(E_theta+E_tc*dthett)+gloc(nphi+i-2,icg)
enddo
C Ufff.... We've done all this!!!
return
enddo
C Ufff.... We've done all this!!!
return
@@
-4862,7
+4940,7
@@
c lprn1=.false.
etheta=etheta+ethetai
if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii
if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1
etheta=etheta+ethetai
if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii
if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1
- gloc(nphi+i-2,icg)=wang*dethetai
+ gloc(nphi+i-2,icg)=wang*dethetai+gloc(nphi+i-2,icg)
enddo
return
end
enddo
return
end
diff --git
a/source/unres/src_MD-M/parmread.F
b/source/unres/src_MD-M/parmread.F
index
9db9b60
..
03efdb2
100644
(file)
--- a/
source/unres/src_MD-M/parmread.F
+++ b/
source/unres/src_MD-M/parmread.F
@@
-991,21
+991,22
@@
c Dtilde(1,1,i)=0.0d0
c Dtilde(1,2,i)=0.0d0
c Dtilde(2,1,i)=0.0d0
c Dtilde(2,2,i)=0.0d0
c Dtilde(1,2,i)=0.0d0
c Dtilde(2,1,i)=0.0d0
c Dtilde(2,2,i)=0.0d0
- EE(1,1,i)= b(10)+b(11)
- EE(2,2,i)=-b(10)+b(11)
- EE(2,1,i)= b(12)-b(13)
- EE(1,2,i)= b(12)+b(13)
- EE(1,1,-i)= b(10)+b(11)
- EE(2,2,-i)=-b(10)+b(11)
- EE(2,1,-i)=-b(12)+b(13)
- EE(1,2,-i)=-b(12)-b(13)
-
+ EEold(1,1,i)= b(10)+b(11)
+ EEold(2,2,i)=-b(10)+b(11)
+ EEold(2,1,i)= b(12)-b(13)
+ EEold(1,2,i)= b(12)+b(13)
+ EEold(1,1,-i)= b(10)+b(11)
+ EEold(2,2,-i)=-b(10)+b(11)
+ EEold(2,1,-i)=-b(12)+b(13)
+ EEold(1,2,-i)=-b(12)-b(13)
+ write(iout,*) "TU DOCHODZE"
c ee(1,1,i)=1.0d0
c ee(2,2,i)=1.0d0
c ee(2,1,i)=0.0d0
c ee(1,2,i)=0.0d0
c ee(2,1,i)=ee(1,2,i)
enddo
c ee(1,1,i)=1.0d0
c ee(2,2,i)=1.0d0
c ee(2,1,i)=0.0d0
c ee(1,2,i)=0.0d0
c ee(2,1,i)=ee(1,2,i)
enddo
+c lprint=.true.
if (lprint) then
do i=1,nloctyp
write (iout,*) 'Type',i
if (lprint) then
do i=1,nloctyp
write (iout,*) 'Type',i
@@
-1023,10
+1024,11
@@
c ee(2,1,i)=ee(1,2,i)
enddo
write(iout,*) 'EE'
do j=1,2
enddo
write(iout,*) 'EE'
do j=1,2
- write (iout,'(2f10.5)') EE(j,1,i),EE(j,2,i)
+ write (iout,'(2f10.5)') EEold(j,1,i),EEold(j,2,i)
enddo
enddo
endif
enddo
enddo
endif
+c lprint=.false.
C
C Read electrostatic-interaction parameters
C
C Read electrostatic-interaction parameters