c call int_from_cart1(.false.)
c endif
#endif
+#ifndef DFA
+ edfadis=0.0d0
+ edfator=0.0d0
+ edfanei=0.0d0
+ edfabet=0.0d0
+#endif
#ifdef TIMING
#ifdef MPI
time00=MPI_Wtime()
C Calculate electrostatic (H-bonding) energy of the main chain.
C
107 continue
+#ifdef DFA
C BARTEK for dfa test!
if (wdfa_dist.gt.0) then
call edfad(edfadis)
else
edfabet=0
endif
+#endif
c print*, 'edfab is finished!', edfabet
cmc
cmc Sep-06: egb takes care of dynamic ss bonds too
if (constr_homology.ge.1) then
call e_modeller(ehomology_constr)
+c print *,'iset=',iset,'me=',me,ehomology_constr,
+c & 'Processor',fg_rank,' CG group',kolor,
+c & ' absolute rank',MyRank
else
ehomology_constr=0.0d0
endif
include 'COMMON.TIME1'
include 'COMMON.MAXGRAD'
include 'COMMON.SCCOR'
+ include 'COMMON.MD'
#ifdef TIMING
#ifdef MPI
time01=MPI_Wtime()
#endif
enddo
enddo
+ if (constr_homology.gt.0) then
+ do i=1,nct
+ do j=1,3
+ gradc(j,i,icg)=gradc(j,i,icg)+duscdiff(j,i)
+ gradx(j,i,icg)=gradx(j,i,icg)+duscdiffx(j,i)
+ enddo
+ enddo
+ endif
#ifdef DEBUG
write (iout,*) "gloc before adding corr"
do i=1,4*nres
& sinph1ph2(maxdouble,maxdouble)
logical lprn /.false./, lprn1 /.false./
etheta=0.0D0
+c write (iout,*) "EBEND ithet_start",ithet_start,
+c & " ithet_end",ithet_end
do i=ithet_start,ithet_end
+ if ((itype(i-1).eq.ntyp1).or.(itype(i-2).eq.ntyp1).or.
+ &(itype(i).eq.ntyp1)) cycle
dethetai=0.0d0
dephii=0.0d0
dephii1=0.0d0
coskt(k)=dcos(k*theti2)
sinkt(k)=dsin(k*theti2)
enddo
- if (i.gt.3) then
+C if (i.gt.3) then
+ if (i.gt.3 .and. itype(i-3).ne.ntyp1) then
#ifdef OSF
phii=phi(i)
if (phii.ne.phii) phii=150.0
enddo
else
phii=0.0d0
- ityp1=nthetyp+1
+ ityp1=ithetyp(itype(i-2))
do k=1,nsingle
cosph1(k)=0.0d0
sinph1(k)=0.0d0
enddo
endif
- if (i.lt.nres) then
+ if ((i.lt.nres).and. itype(i+1).ne.ntyp1) then
#ifdef OSF
phii1=phi(i+1)
if (phii1.ne.phii1) phii1=150.0
enddo
else
phii1=0.0d0
- ityp3=nthetyp+1
+ ityp3=ithetyp(itype(i))
do k=1,nsingle
cosph2(k)=0.0d0
sinph2(k)=0.0d0
& phii1*rad2deg,ethetai
c lprn1=.false.
etheta=etheta+ethetai
+ if (energy_dec) write (iout,'(a6,i5,0pf7.3)')
+ & 'ebend',i,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)=gloc(nphi+i-2,icg)+wang*dethetai
c------------------------------------------------------------------------------
c LICZENIE WIEZOW Z ROWNANIA ENERGII MODELLERA
subroutine e_modeller(ehomology_constr)
- ehomology_constr=0.0
+ ehomology_constr=0.0d0
write (iout,*) "!!!!!UWAGA, JESTEM W DZIWNEJ PETLI, TEST!!!!!"
return
end
dij=dist(i,j)
c write (iout,*) "dij(",i,j,") =",dij
do k=1,constr_homology
+c write(iout,*) ii,k,i,j,l_homo(k,ii),dij,odl(k,ii)
+ if(.not.l_homo(k,ii)) cycle
distance(k)=odl(k,ii)-dij
c write (iout,*) "distance(",k,") =",distance(k)
+c
+c For Gaussian-type Urestr
+c
distancek(k)=0.5d0*distance(k)**2*sigma_odl(k,ii) ! waga_dist rmvd from Gaussian argument
c write (iout,*) "sigma_odl(",k,ii,") =",sigma_odl(k,ii)
c write (iout,*) "distancek(",k,") =",distancek(k)
c distancek(k)=0.5d0*waga_dist*distance(k)**2*sigma_odl(k,ii)
+c
+c For Lorentzian-type Urestr
+c
+ if (waga_dist.lt.0.0d0) then
+ sigma_odlir(k,ii)=dsqrt(1/sigma_odl(k,ii))
+ distancek(k)=distance(k)**2/(sigma_odlir(k,ii)*
+ & (distance(k)**2+sigma_odlir(k,ii)**2))
+ endif
enddo
min_odl=minval(distancek)
odleg2=0.0d0
do k=1,constr_homology
c Nie wiem po co to liczycie jeszcze raz!
-c odleg3=-waga_dist*((distance(i,j,k)**2)/
+c odleg3=-waga_dist(iset)*((distance(i,j,k)**2)/
c & (2*(sigma_odl(i,j,k))**2))
+ if(.not.l_homo(k,ii)) cycle
+ if (waga_dist.ge.0.0d0) then
+c
+c For Gaussian-type Urestr
+c
godl(k)=dexp(-distancek(k)+min_odl)
odleg2=odleg2+godl(k)
+c
+c For Lorentzian-type Urestr
+c
+ else
+ odleg2=odleg2+distancek(k)
+ endif
ccc write(iout,779) i,j,k, "odleg2=",odleg2, "odleg3=", odleg3,
ccc & "dEXP(odleg3)=", dEXP(odleg3),"distance(i,j,k)^2=",
write (iout,*) "godl",(godl(k),k=1,constr_homology) ! exponents
write (iout,*) "ii i j",ii,i,j," odleg2",odleg2 ! sum of exps
#endif
- odleg=odleg-dLOG(odleg2/constr_homology)+min_odl
+ if (waga_dist.ge.0.0d0) then
+c
+c For Gaussian-type Urestr
+c
+ odleg=odleg-dLOG(odleg2/constr_homology)+min_odl
+c
+c For Lorentzian-type Urestr
+c
+ else
+ odleg=odleg+odleg2/constr_homology
+ endif
+c
c write (iout,*) "odleg",odleg ! sum of -ln-s
c Gradient
- sum_godl=odleg2
+c
+c For Gaussian-type Urestr
+c
+ if (waga_dist.ge.0.0d0) sum_godl=odleg2
sum_sgodl=0.0d0
do k=1,constr_homology
c godl=dexp(((-(distance(i,j,k)**2)/(2*(sigma_odl(i,j,k))**2))
c & *waga_dist)+min_odl
c sgodl=-godl(k)*distance(k)*sigma_odl(k,ii)*waga_dist
+c
+ if(.not.l_homo(k,ii)) cycle
+ if (waga_dist.ge.0.0d0) then
+c For Gaussian-type Urestr
+c
sgodl=-godl(k)*distance(k)*sigma_odl(k,ii) ! waga_dist rmvd
+c
+c For Lorentzian-type Urestr
+c
+ else
+ sgodl=-2*sigma_odlir(k,ii)*(distance(k)/(distance(k)**2+
+ & sigma_odlir(k,ii)**2)**2)
+ endif
sum_sgodl=sum_sgodl+sgodl
c sgodl2=sgodl2+sgodl
c write(iout,*) "constr_homology=",constr_homology
c write(iout,*) i, j, k, "TEST K"
enddo
-
- grad_odl3=waga_dist*sum_sgodl/(sum_godl*dij)
+ if (waga_dist.ge.0.0d0) then
+c
+c For Gaussian-type Urestr
+c
+ grad_odl3=waga_homology(iset)*waga_dist
+ & *sum_sgodl/(sum_godl*dij)
+c
+c For Lorentzian-type Urestr
+c
+ else
+c Original grad expr modified by analogy w Gaussian-type Urestr grad
+c grad_odl3=-waga_homology(iset)*waga_dist*sum_sgodl
+ grad_odl3=-waga_homology(iset)*waga_dist*
+ & sum_sgodl/(constr_homology*dij)
+ endif
+c
c grad_odl3=sum_sgodl/(sum_godl*dij)
sum_sgdih=sum_sgdih+sgdih
enddo
c grad_dih3=sum_sgdih/sum_gdih
- grad_dih3=waga_angle*sum_sgdih/sum_gdih
+ grad_dih3=waga_homology(iset)*waga_angle*sum_sgdih/sum_gdih
c write(iout,*)i,k,gdih,sgdih,beta(i+1,i+2,i+3,i+4),grad_dih3
ccc write(iout,747) "GRAD_KAT_1", i, nphi, icg, grad_dih3,
c sgtheta=-gtheta(k)*theta_diff(k)*sigma_theta(k,i)*waga_theta ! right functional form?
sum_sgtheta=sum_sgtheta+sgtheta ! cum variable
enddo
-c grad_theta3=sum_sgtheta/sum_gtheta 1/*theta(i)? s. line below
-c grad_theta3=sum_sgtheta/sum_gtheta
-c
c Final value of gradient using same var as in Econstr_back
- dutheta(i-2)=sum_sgtheta/sum_gtheta*waga_theta
+ gloc(nphi+i-2,icg)=gloc(nphi+i-2,icg)
+ & +sum_sgtheta/sum_gtheta*waga_theta
+ & *waga_homology(iset)
+c dutheta(i-2)=sum_sgtheta/sum_gtheta*waga_theta
+c & *waga_homology(iset)
c dutheta(i)=sum_sgtheta/sum_gtheta
c
c Uconst_back=Uconst_back+waga_theta*utheta(i) ! waga_theta added as weight
c
c
c New implementation
- sum_guscdiff = waga_d*sum_guscdiff
+ sum_guscdiff = waga_homology(iset)*waga_d*sum_guscdiff
do jik=1,3
duscdiff(jik,i-1)=duscdiff(jik,i-1)+
& sum_guscdiff*(dXX_C1tab(jik,i)*dxx+
c Addition of energy of theta angle and SC local geom over constr_homologs ref strs
c
c ehomology_constr=odleg+kat
- ehomology_constr=waga_dist*odleg+waga_angle*kat+waga_theta*Eval
- & +waga_d*Erot
-c write (iout,*) "odleg",odleg," kat",kat," Uconst_back",Uconst_back
-c write (iout,*) "ehomology_constr",ehomology_constr
-c ehomology_constr=odleg+kat+Uconst_back
+c
+c For Lorentzian-type Urestr
+c
+
+ if (waga_dist.ge.0.0d0) then
+c
+c For Gaussian-type Urestr
+c
+ ehomology_constr=(waga_dist*odleg+waga_angle*kat+
+ & waga_theta*Eval+waga_d*Erot)*waga_homology(iset)
+c write (iout,*) "ehomology_constr=",ehomology_constr
+ else
+c
+c For Lorentzian-type Urestr
+c
+ ehomology_constr=(-waga_dist*odleg+waga_angle*kat+
+ & waga_theta*Eval+waga_d*Erot)*waga_homology(iset)
+c write (iout,*) "ehomology_constr=",ehomology_constr
+ endif
+#ifdef DEBUG
+ write (iout,*) "odleg",waga_dist,odleg," kat",waga_angle,kat,
+ & "Eval",waga_theta,eval,
+ & "Erot",waga_d,Erot
+ write (iout,*) "ehomology_constr",ehomology_constr
+#endif
return
c
c FP 01/15 end
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.TORCNSTR'
+ include 'COMMON.CONTROL'
logical lprn
C Set lprn=.true. for debugging
lprn=.false.
c lprn=.true.
etors_d=0.0D0
do i=iphid_start,iphid_end
+ etors_d_ii=0.0D0
itori=itortyp(itype(i-2))
itori1=itortyp(itype(i-1))
itori2=itortyp(itype(i))
sinphi2=dsin(j*phii1)
etors_d=etors_d+v1cij*cosphi1+v1sij*sinphi1+
& v2cij*cosphi2+v2sij*sinphi2
+ if (energy_dec) etors_d_ii=etors_d_ii+
+ & v1cij*cosphi1+v1sij*sinphi1+v2cij*cosphi2+v2sij*sinphi2
gloci1=gloci1+j*(v1sij*cosphi1-v1cij*sinphi1)
gloci2=gloci2+j*(v2sij*cosphi2-v2cij*sinphi2)
enddo
sinphi1m2=dsin(l*phii-(k-l)*phii1)
etors_d=etors_d+v1cdij*cosphi1p2+v2cdij*cosphi1m2+
& v1sdij*sinphi1p2+v2sdij*sinphi1m2
+ if (energy_dec) etors_d_ii=etors_d_ii+
+ & v1cdij*cosphi1p2+v2cdij*cosphi1m2+
+ & v1sdij*sinphi1p2+v2sdij*sinphi1m2
gloci1=gloci1+l*(v1sdij*cosphi1p2+v2sdij*cosphi1m2
& -v1cdij*sinphi1p2-v2cdij*sinphi1m2)
gloci2=gloci2+(k-l)*(v1sdij*cosphi1p2-v2sdij*cosphi1m2
& -v1cdij*sinphi1p2+v2cdij*sinphi1m2)
enddo
enddo
+ if (energy_dec) write (iout,'(a6,i5,0pf7.3)')
+ & 'etor_d',i,etors_d_ii
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)
esccor=0.0D0
do i=itau_start,itau_end
esccor_ii=0.0D0
+ if ((itype(i-2).eq.ntyp1).or.(itype(i-1).eq.ntyp1)) cycle
isccori=isccortyp(itype(i-2))
isccori1=isccortyp(itype(i-1))
phii=phi(i)
projects / unres.git / blobdiff