iatsc_s,iatsc_e,i+1,jj-1,nint_gr(i),istart(i,1),iend(i,1),*12)
ii=nint_gr(i)+1
call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,&
-- iatsc_s,iatsc_e,jj+1,nct,nint_gr(i),istart(i,ii),iend(i,ii),*12)
++ iatsc_s,iatsc_e,jj+1,nct_molec(1),nint_gr(i),istart(i,ii),iend(i,ii),*12)
#else
nint_gr(i)=2
istart(i,1)=i+1
iend(i,1)=jj-1
istart(i,2)=jj+1
-- iend(i,2)=nct
++ iend(i,2)=nct_molec(1)
#endif
endif
else
#ifdef MPI
call int_partition(ind_scint,my_sc_inds,my_sc_inde,i,&
-- iatsc_s,iatsc_e,i+1,nct,nint_gr(i),istart(i,1),iend(i,1),*12)
++ iatsc_s,iatsc_e,i+1,nct_molec(1),nint_gr(i),istart(i,1),iend(i,1),*12)
#else
nint_gr(i)=1
istart(i,1)=i+1
-- iend(i,1)=nct
-- ind_scint=ind_scint+nct-i
++ iend(i,1)=nct_molec(1)
++ ind_scint=ind_scint+nct_molec(1)-i
#endif
endif
#ifdef MPI
iatel_e=0
ind_eleint=0
ind_eleint_old=0
- if (itype(nres_molec(1),1).eq.ntyp_molec(1)) then
- do i=nnt,nct-3
++ if (itype(nres_molec(1),1).eq.ntyp1_molec(1)) then
+ nct_molec(1)=nres_molec(1)-1
+ else
+ nct_molec(1)=nres_molec(1)
+ endif
-
++! print *,"nct",nct,nct_molec(1),itype(nres_molec(1),1),ntyp_molec(1)
+ do i=nnt,nct_molec(1)-3
ijunk=0
call int_partition(ind_eleint,my_ele_inds,my_ele_inde,i,&
-- iatel_s,iatel_e,i+ispp,nct-1,ijunk,ielstart(i),ielend(i),*13)
++ iatel_s,iatel_e,i+ispp,nct_molec(1)-1,ijunk,ielstart(i),ielend(i),*13)
enddo ! i
13 continue
if (iatel_s.eq.0) iatel_s=1
ijunk=0
call int_partition(ind_eleint_vdw,my_ele_inds_vdw,&
my_ele_inde_vdw,i,&
-- iatel_s_vdw,iatel_e_vdw,i+2,nct-1,ijunk,ielstart_vdw(i),&
++ iatel_s_vdw,iatel_e_vdw,i+2,nct_molec(1)-1,ijunk,ielstart_vdw(i),&
ielend_vdw(i),*15)
! write (iout,*) i," ielstart_vdw",ielstart_vdw(i),
! & " ielend_vdw",ielend_vdw(i)
if (i.lt.nnt+iscp) then
!d write (iout,*) 'i.le.nnt+iscp'
call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,&
-- iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,1),&
++ iatscp_s,iatscp_e,i+iscp,nct_molec(1),nscp_gr(i),iscpstart(i,1),&
iscpend(i,1),*14)
else if (i.gt.nct-iscp) then
!d write (iout,*) 'i.gt.nct-iscp'
iscpend(i,1),*14)
ii=nscp_gr(i)+1
call int_partition(ind_scpint,my_scp_inds,my_scp_inde,i,&
-- iatscp_s,iatscp_e,i+iscp,nct,nscp_gr(i),iscpstart(i,ii),&
++ iatscp_s,iatscp_e,i+iscp,nct_molec(1),nscp_gr(i),iscpstart(i,ii),&
iscpend(i,ii),*14)
endif
enddo ! i
iphid_end=iturn4_end+2
iturn4_start=iturn4_start-1
iturn4_end=iturn4_end-1
- print *,"TUTUTU",nres_molec(1),nres
-- call int_bounds(nres-2,ibond_start,ibond_end)
++! print *,"TUTUTU",nres_molec(1),nres
++ call int_bounds(nres_molec(1)-2,ibond_start,ibond_end)
ibond_start=ibond_start+1
ibond_end=ibond_end+1
- call int_bounds(nct-nnt,ibondp_start,ibondp_end)
++ print *,ibond_start,ibond_end
+ call int_bounds(nct_molec(1)-nnt,ibondp_start,ibondp_end)
ibondp_start=ibondp_start+nnt
ibondp_end=ibondp_end+nnt
- call int_bounds1(nres-1,ivec_start,ivec_end)
+ call int_bounds1(nres_molec(1)-1,ivec_start,ivec_end)
! print *,"Processor",myrank,fg_rank,fg_rank1,
! & " ivec_start",ivec_start," ivec_end",ivec_end
iset_start=loc_start+2
! Calculate the bond-stretching energy
!
call ebond(estr)
++ print *,"EBOND",estr
! write(iout,*) "in etotal afer ebond",ipot
!
! endif
enddo
estr=0.5d0*AKP*estr+estr1
++ print *,"estr_bb",estr,AKP
!
! 09/18/07 AL: multimodal bond potential based on AM1 CA-SC PMF's included
!
do i=ibond_start,ibond_end
- iti=iabs(itype(i))
+ iti=iabs(itype(i,1))
++ if (iti.eq.0) print *,"WARNING WRONG SETTTING",i
if (iti.ne.10 .and. iti.ne.ntyp1) then
nbi=nbondterm(iti)
if (nbi.eq.1) then
"estr sc",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
++ print *,"estr_sc",estr
do j=1,3
gradbx(j,i)=AKSC(1,iti)*diff*dc(j,i+nres)/vbld(i+nres)
enddo
usumsqder=usumsqder+ud(j)*uprod2
enddo
estr=estr+uprod/usum
- if (energy_dec) write (iout,*) &
++ print *,"estr_sc",estr,i
++
+ if (energy_dec) write (iout,*) &
"estr sc",i,iti,vbld(i+nres),vbldsc0(1,iti),diff,&
AKSC(1,iti),AKSC(1,iti)*diff*diff
do j=1,3
integer :: i,j,ires,nscat
real(kind=8),dimension(3,20) :: sccor
real(kind=8) :: sccmj
- print *,"I am in sccenter",ires,nscat
++! print *,"I am in sccenter",ires,nscat
do j=1,3
sccmj=0.0D0
do i=1,nscat
ishift=ishift-(ires-ishift+ishift1-ires_old-1)
ires=ires-ishift+ishift1
ires_old=ires
- endif
+ endif
- print *,'atom',ires,atom
++! print *,'atom',ires,atom
if (res.eq.'ACE' .or. res.eq.'NHE') then
- itype(ires)=10
+ itype(ires,1)=10
+ else
+ if (atom.eq.'CA '.or.atom.eq.'N ') then
+ molecule=1
+ itype(ires,molecule)=rescode(ires,res,0,molecule)
+! nres_molec(molecule)=nres_molec(molecule)+1
else
- itype(ires)=rescode(ires,res,0)
+ molecule=2
+ itype(ires,molecule)=rescode(ires,res(3:4),0,molecule)
+! nres_molec(molecule)=nres_molec(molecule)+1
+ endif
endif
else
ires=ires-ishift+ishift1
do j=1,3
sccor(j,iii)=c(j,ires)
enddo
-! write (*,*) card(23:27),ires,itype(ires)
+ else if (.not.unres_pdb .and. (atom.eq."C1'" .or. &
+ atom.eq."C2'" .or. atom.eq."C3'" &
+ .or. atom.eq."C4'" .or. atom.eq."O4'")) then
+ read(card(31:54),'(3f8.3)') (ccc(j),j=1,3)
+!c write (2,'(i5,3f10.5)') ires,(ccc(j),j=1,3)
+ print *,ires,ishift,ishift1
+ counter=counter+1
+! iii=iii+1
+! do j=1,3
+! sccor(j,iii)=c(j,ires)
+! enddo
+ do j=1,3
+ c(j,ires)=c(j,ires)+ccc(j)/5.0
+ enddo
+ if (counter.eq.5) then
+! iii=iii+1
+ nres_molec(molecule)=nres_molec(molecule)+1
+! do j=1,3
+! sccor(j,iii)=c(j,ires)
+! enddo
+ counter=0
+ endif
- print *, "ATOM",atom(1:3)
++! print *, "ATOM",atom(1:3)
+ else if (atom(1:3).eq."C5'") then
+ read (card(19:19),'(a1)') sugar
+ isugar=sugarcode(sugar,ires)
+ if (ibeg.eq.1) then
+ istype(1)=isugar
+ else
+ istype(ires)=isugar
+ endif
+ if (unres_pdb) then
+ read(card(31:54),'(3f8.3)') (c(j,ires),j=1,3)
+ else
+ iii=iii+1
+ read(card(31:54),'(3f8.3)') (sccor(j,iii),j=1,3)
+ endif
+! write (*,*) card(23:27),ires,itype(ires,1)
else if (atom.ne.'O'.and.atom(1:1).ne.'H' .and. &
atom.ne.'N' .and. atom.ne.'C' .and. &
atom(:2).ne.'1H' .and. atom(:2).ne.'2H' .and. &
! nres=ires
nsup=nres
nstart_sup=1
- print *,"molecule",molecule
- if (itype(nres).ne.10) then
++! print *,"molecule",molecule
+ if (itype(nres,1).ne.10) then
nres=nres+1
- itype(nres)=ntyp1
+ itype(nres,molecule)=ntyp1_molec(molecule)
+ nres_molec(molecule)=nres_molec(molecule)+1
if (unres_pdb) then
! 2/15/2013 by Adam: corrected insertion of the last dummy residue
call refsys(nres-3,nres-2,nres-1,e1,e2,e3,fail)
enddo
endif
endif
- print *,'I have',nres, nres_molec(:)
++! print *,'I have',nres, nres_molec(:)
+
!el kontrola nres w pliku inputowym WHAM-a w porownaniu z wartoscia wczytana z pliku pdb
#ifdef WHAM_RUN
if (nres.ne.nres0) then
(c(j,nres+ires),j=1,3)
enddo
endif
+! NOW LETS ROCK! SORTING
+ allocate(c_temporary(3,2*nres))
+ allocate(itype_temporary(nres,5))
+ itype_temporary(:,:)=0
+ seqalingbegin=1
+ do k=1,5
+ do i=1,nres
+ if (itype(i,k).ne.0) then
+ do j=1,3
+ c_temporary(j,seqalingbegin)=c(j,i)
+ c_temporary(j,seqalingbegin+nres)=c(j,i+nres)
+
+ enddo
+ itype_temporary(seqalingbegin,k)=itype(i,k)
+ seqalingbegin=seqalingbegin+1
+ endif
+ enddo
+ enddo
+ do i=1,2*nres
+ do j=1,3
+ c(j,i)=c_temporary(j,i)
+ enddo
+ enddo
+ do k=1,5
+ do i=1,nres
+ itype(i,k)=itype_temporary(i,k)
+ enddo
+ enddo
++ if (itype(1,1).eq.ntyp1) then
++ nsup=nsup-1
++ nstart_sup=2
++ if (unres_pdb) then
++! 2/15/2013 by Adam: corrected insertion of the first dummy residue
++ call refsys(2,3,4,e1,e2,e3,fail)
++ if (fail) then
++ e2(1)=0.0d0
++ e2(2)=1.0d0
++ e2(3)=0.0d0
++ endif
++ do j=1,3
++ c(j,1)=c(j,2)-1.9d0*e2(j)
++ enddo
++ else
++ do j=1,3
++ dcj=(c(j,4)-c(j,3))/2.0
++ c(j,1)=c(j,2)-dcj
++ c(j,nres+1)=c(j,1)
++ enddo
++ endif
++ endif
++
+ if (lprn) then
+ write (iout,'(/a)') &
+ "Cartesian coordinates of the reference structure after sorting"
+ write (iout,'(a,3(3x,a5),5x,3(3x,a5))') &
+ "Residue","X(CA)","Y(CA)","Z(CA)","X(SC)","Y(SC)","Z(SC)"
+ do ires=1,nres
+ write (iout,'(5(a3,1x),i3,3f8.3,5x,3f8.3)') &
+ (restyp(itype(ires,j),j),j=1,5),ires,(c(j,ires),j=1,3),&
+ (c(j,ires+nres),j=1,3)
+ enddo
+ endif
+ print *,seqalingbegin,nres
if(.not.allocated(vbld)) then
allocate(vbld(2*nres))
do i=1,2*nres
projects / unres4.git / commitdiff