X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fcluster%2Fwham%2Fsrc-HCD-5D%2Fenergy_p_new.F;h=5d07d5db79749520c5530904889c4b22b6c7dd5c;hb=79fd686e9d9aaa23e4fda802ee16379470e961d3;hp=5cc851c69766c8c8ab6b4e984967e004da84d3cf;hpb=b6b1679cfe492bd5b73e03410f82917bd9c3aa7b;p=unres.git diff --git a/source/cluster/wham/src-HCD-5D/energy_p_new.F b/source/cluster/wham/src-HCD-5D/energy_p_new.F index 5cc851c..5d07d5d 100644 --- a/source/cluster/wham/src-HCD-5D/energy_p_new.F +++ b/source/cluster/wham/src-HCD-5D/energy_p_new.F @@ -1775,13 +1775,17 @@ C to calculate the el-loc multibody terms of various order. C c write(iout,*) 'SET_MATRICES nphi=',nphi,nres do i=3,nres+1 - if (i.gt. nnt+2 .and. i.lt.nct+2) then + ii=ireschain(i-2) + if (ii.eq.0) cycle + innt=chain_border(1,ii) + inct=chain_border(2,ii) + if (i.gt. innt+2 .and. i.lt.inct+2) then iti = itype2loc(itype(i-2)) else iti=nloctyp endif c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then - if (i.gt. nnt+1 .and. i.lt.nct+1) then + if (i.gt. innt+1 .and. i.lt.inct+1) then iti1 = itype2loc(itype(i-1)) else iti1=nloctyp @@ -1856,6 +1860,19 @@ c b2tilde(2,i-2)=-b2(2,i-2) write (iout,*) 'theta=', theta(i-1) #endif #else + if (i.gt. innt+2 .and. i.lt.inct+2) then +c if (i.gt. nnt+2 .and. i.lt.nct+2) then + iti = itype2loc(itype(i-2)) + else + iti=nloctyp + endif +c write (iout,*) "i",i-1," itype",itype(i-2)," iti",iti +c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then + if (i.gt. nnt+1 .and. i.lt.nct+1) then + iti1 = itype2loc(itype(i-1)) + else + iti1=nloctyp + endif c if (i.gt. nnt+2 .and. i.lt.nct+2) then c iti = itype2loc(itype(i-2)) c else @@ -4520,6 +4537,10 @@ c estr1=0.0d0 c write (iout,*) "distchainmax",distchainmax do i=nnt+1,nct +#ifdef FIVEDIAG + if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) cycle + diff = vbld(i)-vbldp0 +#else if (itype(i-1).eq.ntyp1 .and. itype(i).eq.ntyp1) cycle C estr1=estr1+gnmr1(vbld(i),-1.0d0,distchainmax) C do j=1,3 @@ -4537,6 +4558,9 @@ C write(iout,*) i,diff diff = vbld(i)-vbldp0 c write (iout,*) i,vbld(i),vbldp0,diff,AKP*diff*diff endif +#endif + if (energy_dec) write (iout,'(a7,i5,4f7.3)') + & "estr bb",i,vbld(i),vbldp0,diff,AKP*diff*diff estr=estr+diff*diff do j=1,3 gradb(j,i-1)=AKP*diff*dc(j,i-1)/vbld(i) @@ -4555,8 +4579,9 @@ c nbi=nbondterm(iti) if (nbi.eq.1) then diff=vbld(i+nres)-vbldsc0(1,iti) -C write (iout,*) i,iti,vbld(i+nres),vbldsc0(1,iti),diff, -C & AKSC(1,iti),AKSC(1,iti)*diff*diff + if (energy_dec) write (iout,*) + & i,iti,vbld(i+nres),vbldsc0(1,iti),diff, + & AKSC(1,iti),AKSC(1,iti)*diff*diff estr=estr+0.5d0*AKSC(1,iti)*diff*diff do j=1,3 gradbx(j,i)=AKSC(1,iti)*diff*dc(j,i+nres)/vbld(i+nres) @@ -10178,16 +10203,21 @@ c enddo c min_odl=minval(distancek) - do kk=1,constr_homology - if(l_homo(kk,ii)) then - min_odl=distancek(kk) - exit - endif - enddo - do kk=1,constr_homology - if(l_homo(kk,ii) .and. distancek(kk).lt.min_odl) + if (nexl.gt.0) then + min_odl=0.0d0 + else + do kk=1,constr_homology + if(l_homo(kk,ii)) then + min_odl=distancek(kk) + exit + endif + enddo + do kk=1,constr_homology + if(l_homo(kk,ii) .and. distancek(kk).lt.min_odl) & min_odl=distancek(kk) - enddo + enddo + endif + c write (iout,* )"min_odl",min_odl #ifdef DEBUG write (iout,*) "ij dij",i,j,dij