X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD%2Fenergy_p_new_barrier.F;h=4bccdf36ae8d1a6d9fdb3585da5079c52a4ffc18;hb=34d3ad3987785642be58fb2f26557d3314215577;hp=e3c154ea932fe54e6ea5621568cbd742a7b5d947;hpb=f754eecac3c8a4547d905cca354f9a0db1531729;p=unres.git diff --git a/source/unres/src_MD/energy_p_new_barrier.F b/source/unres/src_MD/energy_p_new_barrier.F index e3c154e..4bccdf3 100644 --- a/source/unres/src_MD/energy_p_new_barrier.F +++ b/source/unres/src_MD/energy_p_new_barrier.F @@ -131,6 +131,11 @@ C C Calculate electrostatic (H-bonding) energy of the main chain. C 107 continue +cmc +cmc Sep-06: egb takes care of dynamic ss bonds too +cmc +c if (dyn_ss) call dyn_set_nss + c print *,"Processor",myrank," computed USCSC" #ifdef TIMING #ifdef MPI @@ -326,6 +331,7 @@ C energia(23)=evdw_m c print *," Processor",myrank," calls SUM_ENERGY" call sum_energy(energia,.true.) + if (dyn_ss) call dyn_set_nss c print *," Processor",myrank," left SUM_ENERGY" #ifdef TIMING #ifdef MPI @@ -423,14 +429,14 @@ cMS$ATTRIBUTES C :: proc_proc #ifdef SPLITELE etot=wsc*evdw+wscp*evdw2+welec*ees+wvdwpp*evdw1 & +wang*ebe+wtor*etors+wscloc*escloc - & +wstrain*ehpb+nss*ebr+wcorr*ecorr+wcorr5*ecorr5 + & +wstrain*ehpb+wcorr*ecorr+wcorr5*ecorr5 & +wcorr6*ecorr6+wturn4*eello_turn4+wturn3*eello_turn3 & +wturn6*eturn6+wel_loc*eel_loc+edihcnstr+wtor_d*etors_d & +wbond*estr+Uconst+wsccor*esccor #else etot=wsc*evdw+wscp*evdw2+welec*(ees+evdw1) & +wang*ebe+wtor*etors+wscloc*escloc - & +wstrain*ehpb+nss*ebr+wcorr*ecorr+wcorr5*ecorr5 + & +wstrain*ehpb+wcorr*ecorr+wcorr5*ecorr5 & +wcorr6*ecorr6+wturn4*eello_turn4+wturn3*eello_turn3 & +wturn6*eturn6+wel_loc*eel_loc+edihcnstr+wtor_d*etors_d & +wbond*estr+Uconst+wsccor*esccor @@ -471,7 +477,7 @@ cMS$ATTRIBUTES C :: proc_proc 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) include 'COMMON.SETUP' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' @@ -483,6 +489,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 +762,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" @@ -774,6 +780,19 @@ c enddo do i=1,4*nres glocbuf(i)=gloc(i,icg) enddo +#ifdef DEBUG + write (iout,*) "gloc_sc before reduce" + do i=1,nres + do j=1,3 + write (iout,*) i,j,gloc_sc(j,i,icg) + enddo + enddo +#endif + 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 @@ -784,8 +803,18 @@ c enddo & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR) call MPI_Reduce(glocbuf(1),gloc(1,icg),4*nres, & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR) + 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 #ifdef DEBUG + write (iout,*) "gloc_sc after reduce" + do i=1,nres + do j=1,3 + write (iout,*) i,j,gloc_sc(j,i,icg) + enddo + enddo +#endif +#ifdef DEBUG write (iout,*) "gloc after reduce" do i=1,4*nres write (iout,*) i,gloc(i,icg) @@ -1029,25 +1058,25 @@ C------------------------------------------------------------------------ & edihcnstr,ebr*nss, & Uconst,etot 10 format (/'Virtual-chain energies:'// - & 'EVDW= ',1pE16.6,' WEIGHT=',1pD16.6,' (SC-SC)'/ - & 'EVDW2= ',1pE16.6,' WEIGHT=',1pD16.6,' (SC-p)'/ - & 'EES= ',1pE16.6,' WEIGHT=',1pD16.6,' (p-p)'/ - & 'EVDWPP=',1pE16.6,' WEIGHT=',1pD16.6,' (p-p VDW)'/ - & 'ESTR= ',1pE16.6,' WEIGHT=',1pD16.6,' (stretching)'/ - & 'EBE= ',1pE16.6,' WEIGHT=',1pD16.6,' (bending)'/ - & 'ESC= ',1pE16.6,' WEIGHT=',1pD16.6,' (SC local)'/ - & 'ETORS= ',1pE16.6,' WEIGHT=',1pD16.6,' (torsional)'/ - & 'ETORSD=',1pE16.6,' WEIGHT=',1pD16.6,' (double torsional)'/ - & 'EHBP= ',1pE16.6,' WEIGHT=',1pD16.6, + & 'EVDW= ',1pE16.6,' WEIGHT=',1pE16.6,' (SC-SC)'/ + & 'EVDW2= ',1pE16.6,' WEIGHT=',1pE16.6,' (SC-p)'/ + & 'EES= ',1pE16.6,' WEIGHT=',1pE16.6,' (p-p)'/ + & 'EVDWPP=',1pE16.6,' WEIGHT=',1pE16.6,' (p-p VDW)'/ + & 'ESTR= ',1pE16.6,' WEIGHT=',1pE16.6,' (stretching)'/ + & 'EBE= ',1pE16.6,' WEIGHT=',1pE16.6,' (bending)'/ + & 'ESC= ',1pE16.6,' WEIGHT=',1pE16.6,' (SC local)'/ + & 'ETORS= ',1pE16.6,' WEIGHT=',1pE16.6,' (torsional)'/ + & 'ETORSD=',1pE16.6,' WEIGHT=',1pE16.6,' (double torsional)'/ + & 'EHPB= ',1pE16.6,' WEIGHT=',1pE16.6, & ' (SS bridges & dist. cnstr.)'/ - & 'ECORR4=',1pE16.6,' WEIGHT=',1pD16.6,' (multi-body)'/ - & 'ECORR5=',1pE16.6,' WEIGHT=',1pD16.6,' (multi-body)'/ - & 'ECORR6=',1pE16.6,' WEIGHT=',1pD16.6,' (multi-body)'/ - & 'EELLO= ',1pE16.6,' WEIGHT=',1pD16.6,' (electrostatic-local)'/ - & 'ETURN3=',1pE16.6,' WEIGHT=',1pD16.6,' (turns, 3rd order)'/ - & 'ETURN4=',1pE16.6,' WEIGHT=',1pD16.6,' (turns, 4th order)'/ - & 'ETURN6=',1pE16.6,' WEIGHT=',1pD16.6,' (turns, 6th order)'/ - & 'ESCCOR=',1pE16.6,' WEIGHT=',1pD16.6,' (backbone-rotamer corr)'/ + & 'ECORR4=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ + & 'ECORR5=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ + & 'ECORR6=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/ + & 'EELLO= ',1pE16.6,' WEIGHT=',1pE16.6,' (electrostatic-local)'/ + & 'ETURN3=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 3rd order)'/ + & 'ETURN4=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 4th order)'/ + & 'ETURN6=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 6th order)'/ + & 'ESCCOR=',1pE16.6,' WEIGHT=',1pE16.6,' (backbone-rotamer corr)'/ & 'EDIHC= ',1pE16.6,' (dihedral angle constraints)'/ & 'ESS= ',1pE16.6,' (disulfide-bridge intrinsic energy)'/ & 'UCONST= ',1pE16.6,' (Constraint energy)'/ @@ -1529,6 +1558,7 @@ C include 'COMMON.IOUNITS' include 'COMMON.CALC' include 'COMMON.CONTROL' + include 'COMMON.SBRIDGE' logical lprn evdw=0.0D0 ccccc energy_dec=.false. @@ -1557,6 +1587,12 @@ C Calculate SC interaction energy. C do iint=1,nint_gr(i) do j=istart(i,iint),iend(i,iint) + IF (dyn_ss_mask(i).and.dyn_ss_mask(j)) THEN + call dyn_ssbond_ene(i,j,evdwij) + evdw=evdw+evdwij + if (energy_dec) write (iout,'(a6,2i5,0pf7.3,a3)') + & 'evdw',i,j,evdwij,' ss' + ELSE ind=ind+1 itypj=itype(j) c dscj_inv=dsc_inv(itypj) @@ -1643,9 +1679,10 @@ c & " eps3rt",eps3rt," eps1",eps1," e1",e1," e2",e2 & evdwij endif - if (energy_dec) write (iout,'(a6,2i5,0pf7.3)') - & 'evdw',i,j,evdwij - + if (energy_dec) then + write (iout,'(a6,2i5,0pf7.3)') 'evdw',i,j,evdwij + call flush(iout) + endif C Calculate gradient components. e1=e1*eps1*eps2rt**2*eps3rt**2 fac=-expon*(e1+evdwij)*rij_shift @@ -1666,6 +1703,7 @@ C Calculate angular part of the gradient. #else call sc_grad #endif + ENDIF ! dyn_ss enddo ! j enddo ! iint enddo ! i @@ -2986,6 +3024,9 @@ C C Loop over i,i+2 and i,i+3 pairs of the peptide groups C do i=iturn3_start,iturn3_end +C if (itype(i).eq.21 .or. itype(i+1).eq.21 +C & .or. itype(i+2).eq.21 .or. itype(i+3).eq.21.or.itype(i+4).eq.21) +C & cycle dxi=dc(1,i) dyi=dc(2,i) dzi=dc(3,i) @@ -3001,6 +3042,10 @@ C num_cont_hb(i)=num_conti enddo do i=iturn4_start,iturn4_end +C if (itype(i).eq.21 .or. itype(i+1).eq.21 +C & .or. itype(i+2).eq.21 .or. itype(i+3).eq.21.or.itype(i+4).eq.21 +C & .or. itype(i+5).eq.21) +C & cycle dxi=dc(1,i) dyi=dc(2,i) dzi=dc(3,i) @@ -3019,6 +3064,8 @@ c c Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3 c do i=iatel_s,iatel_e +C if (itype(i).eq.21 .or. itype(i+1).eq.21 +C &.or.itype(i+2)) cycle dxi=dc(1,i) dyi=dc(2,i) dzi=dc(3,i) @@ -3031,6 +3078,8 @@ c c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i) num_conti=num_cont_hb(i) do j=ielstart(i),ielend(i) +C if (itype(j).eq.21 .or. itype(j+1).eq.21 +C &.or.itype(j+2)) cycle call eelecij(i,j,ees,evdw1,eel_loc) enddo ! j num_cont_hb(i)=num_conti @@ -4243,9 +4292,15 @@ c write (iout,*) "i",i," ii",ii," iii",iii," jj",jj," jjj",jjj, c & dhpb(i),dhpb1(i),forcon(i) C 24/11/03 AL: SS bridges handled separately because of introducing a specific C distance and angle dependent SS bond potential. - if (ii.gt.nres .and. itype(iii).eq.1 .and. itype(jjj).eq.1) then +cmc if (ii.gt.nres .and. itype(iii).eq.1 .and. itype(jjj).eq.1) then +C 18/07/06 MC: Use the convention that the first nss pairs are SS bonds + if (.not.dyn_ss .and. i.le.nss) then +C 15/02/13 CC dynamic SSbond - additional check + if (ii.gt.nres + & .and. itype(iii).eq.1 .and. itype(jjj).eq.1) then call ssbond_ene(iii,jjj,eij) ehpb=ehpb+2*eij + endif cd write (iout,*) "eij",eij else if (ii.gt.nres .and. jj.gt.nres) then c Restraints from contact prediction @@ -4384,7 +4439,7 @@ c dscj_inv=dsc_inv(itypj) deltat12=om2-om1+2.0d0 cosphi=om12-om1*om2 eij=akcm*deltad*deltad+akth*(deltat1*deltat1+deltat2*deltat2) - & +akct*deltad*deltat12 + & +akct*deltad*deltat12+ebr & +v1ss*cosphi+v2ss*cosphi*cosphi+v3ss*cosphi*cosphi*cosphi c write(iout,*) i,j,"rij",rij,"d0cm",d0cm," akcm",akcm," akth",akth, c & " akct",akct," deltad",deltad," deltat",deltat1,deltat2, @@ -4439,6 +4494,8 @@ c do i=ibondp_start,ibondp_end diff = vbld(i)-vbldp0 c write (iout,*) i,vbld(i),vbldp0,diff,AKP*diff*diff + 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) @@ -4457,6 +4514,12 @@ c 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) then + write (iout,*) + & "estr sc",i,iti,vbld(i+nres),vbldsc0(1,iti),diff, + & AKSC(1,iti),AKSC(1,iti)*diff*diff + call flush(iout) + endif 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) @@ -4735,6 +4798,8 @@ C logical lprn /.false./, lprn1 /.false./ etheta=0.0D0 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 @@ -4744,7 +4809,8 @@ C 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(max0(i-3,1)).ne.ntyp1) then #ifdef OSF phii=phi(i) if (phii.ne.phii) phii=150.0 @@ -4758,13 +4824,13 @@ C 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 @@ -4779,7 +4845,7 @@ C enddo else phii1=0.0d0 - ityp3=nthetyp+1 + ityp3=ithetyp(itype(i)) do k=1,nsingle cosph2(k)=0.0d0 sinph2(k)=0.0d0 @@ -4885,10 +4951,12 @@ C enddo enddo 10 continue - if (lprn1) write (iout,'(i2,3f8.1,9h ethetai ,f10.5)') - & i,theta(i)*rad2deg,phii*rad2deg, + if (lprn1) write (iout,'(a4,i2,3f8.1,9h ethetai ,f10.5)') + & 'ebe', i,theta(i)*rad2deg,phii*rad2deg, & phii1*rad2deg,ethetai 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)=wang*dethetai @@ -5691,7 +5759,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 @@ -5806,6 +5874,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) @@ -5829,12 +5898,14 @@ C 6/23/01 Compute double torsional energy 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)) @@ -5842,7 +5913,6 @@ c lprn=.true. phii1=phi(i+1) gloci1=0.0D0 gloci2=0.0D0 -C Regular cosine and sine terms do j=1,ntermd_1(itori,itori1,itori2) v1cij=v1c(1,j,itori,itori1,itori2) v1sij=v1s(1,j,itori,itori1,itori2) @@ -5854,6 +5924,8 @@ C Regular cosine and sine terms 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 @@ -5869,14 +5941,20 @@ C Regular cosine and sine terms 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) enddo return end @@ -5909,26 +5987,66 @@ C Set lprn=.true. for debugging c lprn=.true. c write (iout,*) "EBACK_SC_COR",iphi_start,iphi_end,nterm_sccor esccor=0.0D0 - do i=iphi_start,iphi_end - esccor_ii=0.0D0 - itori=itype(i-2) - itori1=itype(i-1) + do i=itau_start,itau_end +C do i=42,42 + + 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) + +cccc Added 9 May 2012 +cc Tauangle is torsional engle depending on the value of first digit +c(see comment below) +cc Omicron is flat angle depending on the value of first digit +c(see comment below) +C print *,i,tauangle(1,i) + + do intertyp=1,3 !intertyp + esccor_ii=0.0D0 +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.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.eq.4).and.(itype(1).eq.21)) cycle + if ((intertyp.eq.1).and.(i.eq.nres).and.(itype(nres).eq.21)) + & 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)) + if (energy_dec) esccor_ii=esccor_ii+v1ij*cosphi+v2ij*sinphi esccor=esccor+v1ij*cosphi+v2ij*sinphi gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi) enddo + if (energy_dec) write (iout,'(a6,i5,i2,0pf7.3)') + & 'esccor',i,intertyp,esccor_ii +C print *,i,tauangle(1,i),gloci + gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci +c write (iout,*) "WTF",intertyp,i,itype(i),v1ij*cosphi+v2ij*sinphi +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,itori,itori1),j=1,6),(v2sccor(j,itori,itori1),j=1,6) + & (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 +c do i=1,nres +c write (iout,*) "W@T@F", gloc_sc(1,i,icg),gloc_sc(2,i,icg), +c & gloc_sc(3,i,icg) +c enddo return end c---------------------------------------------------------------------------- @@ -8035,9 +8153,9 @@ C C Parallel Antiparallel C C o o -C /l\ /j\ -C / \ / \ -C /| o | | o |\ +C /l\ /j\ +C / \ / \ +C /| o | | o |\ C \ j|/k\| / \ |/k\|l / C \ / \ / \ / \ / C o o o o @@ -8132,22 +8250,22 @@ c---------------------------------------------------------------------------- include 'COMMON.GEO' logical swap double precision vv(2),pizda(2,2),auxmat(2,2),auxvec(2), - & auxvec1(2),auxvec2(1),auxmat1(2,2) + & auxvec1(2),auxvec2(2),auxmat1(2,2) logical lprn common /kutas/ lprn CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C \ /l\ /j\ / -C \ / \ / \ / -C o| o | | o |o -C \ j|/k\| \ |/k\|l -C \ / \ \ / \ -C o o -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C \ /l\ /j\ / C +C \ / \ / \ / C +C o| o | | o |o C +C \ j|/k\| \ |/k\|l C +C \ / \ \ / \ C +C o o C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC cd write (2,*) 'eello6_graph2: i,',i,' j',j,' k',k,' l',l C AL 7/4/01 s1 would occur in the sixth-order moment, @@ -8318,18 +8436,18 @@ c---------------------------------------------------------------------------- double precision vv(2),pizda(2,2),auxmat(2,2),auxvec(2) logical swap CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C /l\ / \ /j\ -C / \ / \ / \ -C /| o |o o| o |\ -C j|/k\| / |/k\|l / -C / \ / / \ / -C / o / o -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C /l\ / \ /j\ C +C / \ / \ / \ C +C /| o |o o| o |\ C +C j|/k\| / |/k\|l / C +C / \ / / \ / C +C / o / o C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C 4/7/01 AL Component s1 was removed, because it pertains to the respective @@ -8435,18 +8553,18 @@ c---------------------------------------------------------------------------- & auxvec1(2),auxmat1(2,2) logical swap CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC -C -C Parallel Antiparallel -C -C o o -C /l\ / \ /j\ -C / \ / \ / \ -C /| o |o o| o |\ -C \ j|/k\| \ |/k\|l -C \ / \ \ / \ -C o \ o \ -C i i -C +C C +C Parallel Antiparallel C +C C +C o o C +C /l\ / \ /j\ C +C / \ / \ / \ C +C /| o |o o| o |\ C +C \ j|/k\| \ |/k\|l C +C \ / \ \ / \ C +C o \ o \ C +C i i C +C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C 4/7/01 AL Component s1 was removed, because it pertains to the respective