--- /dev/null
+ subroutine enecalc(islice,*)
+ implicit none
+ include "DIMENSIONS"
+ include "DIMENSIONS.ZSCOPT"
+ include "DIMENSIONS.FREE"
+#ifdef MPI
+ include "mpif.h"
+ integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+ include "COMMON.MPI"
+#endif
+ include "COMMON.CHAIN"
+ include "COMMON.IOUNITS"
+ include "COMMON.PROTFILES"
+ include "COMMON.NAMES"
+ include "COMMON.VAR"
+ include "COMMON.SBRIDGE"
+ include "COMMON.GEO"
+ include "COMMON.FFIELD"
+ include "COMMON.ENEPS"
+ include "COMMON.LOCAL"
+ include "COMMON.WEIGHTS"
+ include "COMMON.INTERACT"
+ include "COMMON.FREE"
+ include "COMMON.ENERGIES"
+ include "COMMON.CONTROL"
+ include "COMMON.TORCNSTR"
+ character*64 nazwa
+ character*80 bxname
+ character*3 liczba
+ double precision qwolynes
+ external qwolynes
+ integer errmsg_count,maxerrmsg_count /100/
+ double precision rmsnat,gyrate
+ external rmsnat,gyrate
+ double precision tole /1.0d-1/
+ integer i,itj,ii,iii,j,k,l,licz
+ integer ir,ib,ipar,iparm
+ integer iscor,islice
+ real*4 csingle(3,maxres2)
+ double precision energ
+ integer ilen,iroof
+ external ilen,iroof
+ double precision energia(0:max_ene),rmsdev,efree,eini
+ double precision fT(6),quot,quotl,kfacl,kfac /2.4d0/,T0 /3.0d2/
+ double precision tt
+ integer snk_p(MaxR,MaxT_h,Max_parm)
+ logical lerr
+ character*64 bprotfile_temp
+ call opentmp(islice,ientout,bprotfile_temp)
+ iii=0
+ ii=0
+ errmsg_count=0
+ write (iout,*) "enecalc: nparmset ",nparmset
+#ifdef MPI
+ do iparm=1,nParmSet
+ do ib=1,nT_h(iparm)
+ do i=1,nR(ib,iparm)
+ snk_p(i,ib,iparm)=0
+ enddo
+ enddo
+ enddo
+ do i=indstart(me1),indend(me1)
+#else
+ do iparm=1,nParmSet
+ do ib=1,nT_h(iparm)
+ do i=1,nR(ib,iparm)
+ snk(i,ib,iparm)=0
+ enddo
+ enddo
+ enddo
+ do i=1,ntot
+#endif
+ read(ientout,rec=i,err=101)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+ & eini,efree,rmsdev,(q(j,iii+1),j=1,nQ),iR,ib,ipar
+ if (indpdb.gt.0) then
+ do k=1,nres
+ do l=1,3
+ c(l,k)=csingle(l,k)
+ enddo
+ enddo
+ do k=nnt,nct
+ do l=1,3
+ c(l,k+nres)=csingle(l,k+nres)
+ enddo
+ enddo
+ q(nQ+1,iii+1)=rmsnat(iii+1)
+ endif
+ q(nQ+2,iii+1)=gyrate(iii+1)
+c fT=T0*beta_h(ib,ipar)*1.987D-3
+c ft=2.0d0/(1.0d0+1.0d0/(T0*beta_h(ib,ipar)*1.987D-3))
+ if (rescale_mode.eq.1) then
+ quot=1.0d0/(T0*beta_h(ib,ipar)*1.987D-3)
+#if defined(FUNCTH)
+ tt = 1.0d0/(beta_h(ib,ipar)*1.987D-3)
+ ft(6)=(320.0+80.0*dtanh((tt-320.0)/80.0))/320.0
+#elif defined(FUNCT)
+ ft(6)=quot
+#else
+ ft(6)=1.0d0
+#endif
+ quotl=1.0d0
+ kfacl=1.0d0
+ do l=1,5
+ quotl=quotl*quot
+ kfacl=kfacl*kfac
+ fT(l)=kfacl/(kfacl-1.0d0+quotl)
+ enddo
+ else if (rescale_mode.eq.2) then
+ quot=1.0d0/(T0*beta_h(ib,ipar)*1.987D-3)
+#if defined(FUNCTH)
+ tt = 1.0d0/(beta_h(ib,ipar)*1.987D-3)
+ ft(6)=(320.0+80.0*dtanh((tt-320.0)/80.0))/320.0
+#elif defined(FUNCT)
+ ft(6)=quot
+#else
+ ft(6)=1.0d0
+#endif
+ quotl=1.0d0
+ do l=1,5
+ quotl=quotl*quot
+ fT(l)=1.12692801104297249644d0/
+ & dlog(dexp(quotl)+dexp(-quotl))
+ enddo
+ else if (rescale_mode.eq.0) then
+ do l=1,5
+ fT(l)=1.0d0
+ enddo
+ else
+ write (iout,*) "Error in ECECALC: wrong RESCALE_MODE",
+ & rescale_mode
+ call flush(iout)
+ return1
+ endif
+
+c write (iout,*) "T",1.0d0/(beta_h(ib,ipar)*1.987D-3)," T0",T0,
+c & " kfac",kfac,"quot",quot," fT",fT
+ do j=1,2*nres
+ do k=1,3
+ c(k,j)=csingle(k,j)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ ii=ii+1
+ do iparm=1,nparmset
+
+ call restore_parm(iparm)
+#ifdef DEBUG
+ write (iout,*) wsc,wscp,welec,wvdwpp,wang,wtor,wscloc,
+ & wcorr,wcorr5,wcorr6,wturn4,wturn3,wturn6,wel_loc,
+ & wtor_d,wsccor,wbond
+#endif
+c write (iout,*) "Calling ETOTAL"
+ call etotal(energia(0),fT,beta_h(ib,iparm))
+#ifdef DEBUG
+ write (iout,*) "Conformation",i
+ call enerprint(energia(0),fT)
+c write (iout,'(2i5,21f8.2)') i,iparm,(energia(k),k=1,21)
+c write (iout,*) "ftors",ftors
+c call intout
+#endif
+ if (energia(0).ge.1.0d20) then
+ write (iout,*) "NaNs detected in some of the energy",
+ & " components for conformation",ii+1
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+c call intout
+c call pdbout(ii+1,beta_h(ib,ipar),efree,energia(0),0.0d0,rmsdev)
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,*) "The components of the energy are:"
+ call enerprint(energia(0),fT)
+ write (iout,*)
+ & "This conformation WILL NOT be added to the database."
+ call flush(iout)
+ goto 121
+ else
+#ifdef DEBUG
+ if (ipar.eq.iparm) write (iout,*) i,iparm,
+ & 1.0d0/(beta_h(ib,ipar)*1.987D-3),eini,energia(0)
+#endif
+ if (ipar.eq.iparm .and. einicheck.gt.0 .and.
+ & dabs(eini-energia(0)).gt.tole) then
+ if (errmsg_count.le.maxerrmsg_count) then
+ write (iout,'(2a,2e15.5,a,2i8,a,f8.1)')
+ & "Warning: energy differs remarkably from ",
+ & " the value read in: ",energia(0),eini," point",
+ & iii+1,indstart(me1)+iii," T",
+ & 1.0d0/(1.987D-3*beta_h(ib,ipar))
+ errmsg_count=errmsg_count+1
+ if (errmsg_count.gt.maxerrmsg_count)
+ & write (iout,*) "Too many warning messages"
+ if (einicheck.gt.1) then
+ write (iout,*) "Calculation stopped."
+ call flush(iout)
+#ifdef MPI
+ call MPI_Abort(WHAM_COMM,IERROR,ERRCODE)
+#endif
+ call flush(iout)
+ return1
+ endif
+ endif
+ endif
+ potE(iii+1,iparm)=energia(0)
+ do k=1,21
+ enetb(k,iii+1,iparm)=energia(k)
+ enddo
+#ifdef DEBUG
+ write (iout,'(2i5,f10.1,3e15.5)') i,iii,
+ & 1.0d0/(beta_h(ib,ipar)*1.987D-3),energia(0),eini,efree
+ call enerprint(energia(0),fT)
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+ write (iout,'(8f10.5)') (q(k,iii+1),k=1,nQ)
+ write (iout,'(f10.5,i10)') rmsdev,iscor
+ call enerprint(energia(0),fT)
+ write(liczba,'(bz,i3.3)') me
+ nazwa="test"//liczba//".pdb"
+ write (iout,*) "pdb file",nazwa
+ open (ipdb,file=nazwa,position="append")
+ call pdbout(ii+1,beta_h(ib,ipar),efree,energia(0),0.0d0,rmsdev)
+ close(ipdb)
+#endif
+ endif
+
+ enddo ! iparm
+
+ iii=iii+1
+ if (q(1,iii).le.0.0d0 .and. indpdb.gt.0) q(1,iii)=qwolynes(0,0)
+ write (ientout,rec=iii)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+ & potE(iii,ipar),efree,rmsdev,(q(k,iii),k=1,nQ),iR,ib,ipar
+c write (iout,'(2i5,2e15.5)') ii,iii,potE(iii,ipar),efree
+#ifdef MPI
+ if (separate_parset) then
+ snk_p(iR,ib,1)=snk_p(iR,ib,1)+1
+ else
+ snk_p(iR,ib,ipar)=snk_p(iR,ib,ipar)+1
+ endif
+c write (iout,*) "iii",iii," iR",iR," ib",ib," ipar",ipar,
+c & " snk",snk_p(iR,ib,ipar)
+#else
+ snk(iR,ib,ipar,islice)=snk(iR,ib,ipar,islice)+1
+#endif
+ 121 continue
+ enddo
+#ifdef MPI
+ scount(me)=iii
+ write (iout,*) "Me",me," scount",scount(me)
+ call flush(iout)
+c Master gathers updated numbers of conformations written by all procs.
+ call MPI_AllGather( scount(me), 1, MPI_INTEGER, scount(0), 1,
+ & MPI_INTEGER, WHAM_COMM, IERROR)
+ indstart(0)=1
+ indend(0)=scount(0)
+ do i=1, Nprocs-1
+ indstart(i)=indend(i-1)+1
+ indend(i)=indstart(i)+scount(i)-1
+ enddo
+ write (iout,*)
+ write (iout,*) "Revised conformation counts"
+ do i=0,nprocs1-1
+ write (iout,'(a,i5,a,i7,a,i7,a,i7)')
+ & "Processor",i," indstart",indstart(i),
+ & " indend",indend(i)," count",scount(i)
+ enddo
+ call flush(iout)
+ call MPI_AllReduce(snk_p(1,1,1),snk(1,1,1,islice),
+ & MaxR*MaxT_h*nParmSet,
+ & MPI_INTEGER,MPI_SUM,WHAM_COMM,IERROR)
+#endif
+ stot(islice)=0
+ do iparm=1,nParmSet
+ do ib=1,nT_h(iparm)
+ do i=1,nR(ib,iparm)
+ stot(islice)=stot(islice)+snk(i,ib,iparm,islice)
+ enddo
+ enddo
+ enddo
+ write (iout,*) "Revised SNK"
+ do iparm=1,nParmSet
+ do ib=1,nT_h(iparm)
+ write (iout,'("Param",i3," Temp",f6.1,3x,32i8)')
+ & iparm,1.0d0/(1.987D-3*beta_h(ib,iparm)),
+ & (snk(i,ib,iparm,islice),i=1,nR(ib,iparm))
+ write (iout,*) "snk_p",(snk_p(i,ib,iparm),i=1,nR(ib,iparm))
+ enddo
+ enddo
+ write (iout,'("Total",i10)') stot(islice)
+ call flush(iout)
+ return
+ 101 write (iout,*) "Error in scratchfile."
+ call flush(iout)
+ return1
+ end
+c------------------------------------------------------------------------------
+ subroutine write_dbase(islice,*)
+ implicit none
+ include "DIMENSIONS"
+ include "DIMENSIONS.ZSCOPT"
+ include "DIMENSIONS.FREE"
+ include "DIMENSIONS.COMPAR"
+#ifdef MPI
+ include "mpif.h"
+ integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+ include "COMMON.MPI"
+#endif
+ include "COMMON.CONTROL"
+ include "COMMON.CHAIN"
+ include "COMMON.IOUNITS"
+ include "COMMON.PROTFILES"
+ include "COMMON.NAMES"
+ include "COMMON.VAR"
+ include "COMMON.SBRIDGE"
+ include "COMMON.GEO"
+ include "COMMON.FFIELD"
+ include "COMMON.ENEPS"
+ include "COMMON.LOCAL"
+ include "COMMON.WEIGHTS"
+ include "COMMON.INTERACT"
+ include "COMMON.FREE"
+ include "COMMON.ENERGIES"
+ include "COMMON.COMPAR"
+ include "COMMON.PROT"
+ character*64 nazwa
+ character*80 bxname,cxname
+ character*64 bprotfile_temp
+ character*3 liczba,licz
+ character*2 licz2
+ integer i,itj,ii,iii,j,k,l
+ integer ixdrf,iret
+ integer iscor,islice
+ double precision rmsdev,efree,eini
+ real*4 csingle(3,maxres2)
+ double precision energ
+ integer ilen,iroof
+ external ilen,iroof
+ integer ir,ib,iparm
+ write (licz2,'(bz,i2.2)') islice
+ call opentmp(islice,ientout,bprotfile_temp)
+ write (iout,*) "bprotfile_temp ",bprotfile_temp
+ call flush(iout)
+ if (.not.bxfile .and. .not. cxfile .and. indpdb.eq.0
+ & .and. ensembles.eq.0) then
+ close(ientout,status="delete")
+ return
+ endif
+#ifdef MPI
+ write (liczba,'(bz,i3.3)') me
+ if (bxfile .or. cxfile .or. ensembles.gt.0) then
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//liczba//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",
+ & form="unformatted",access="direct",recl=lenrec1)
+ endif
+#else
+ if (bxfile .or. cxfile .or. ensembles.gt.0) then
+ if (nslice.eq.1) then
+ bxname = prefix(:ilen(prefix))//".bx"
+ else
+ bxname = prefix(:ilen(prefix))//
+ & "_slice_"//licz2//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",
+ & form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Calculating energies; writing geometry",
+ & " and energy components to ",bxname(:ilen(bxname))
+ endif
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfopen_(ixdrf,cxname, "w", iret)
+#else
+ call xdrfopen(ixdrf,cxname, "w", iret)
+#endif
+ if (iret.eq.0) then
+ write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
+ cxfile=.fale.
+ endif
+ endif
+#endif
+ if (indpdb.gt.0) then
+ if (nslice.eq.1) then
+#ifdef MPI
+ if (.not.separate_parset) then
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))
+ & //liczba//'.stat'
+ else
+ write (licz,'(bz,i3.3)') myparm
+ statname=prefix(:ilen(prefix))//'_par'//licz//'_'//
+ & pot(:ilen(pot))//liczba//'.stat'
+ endif
+
+#else
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//'.stat'
+#endif
+ else
+#ifdef MPI
+ if (.not.separate_parset) then
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//
+ & "_slice_"//licz2//liczba//'.stat'
+ else
+ write (licz,'(bz,i3.3)') myparm
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))//
+ & '_par'//licz//"_slice_"//licz2//liczba//'.stat'
+ endif
+#else
+ statname=prefix(:ilen(prefix))//'_'//pot(:ilen(pot))
+ & //"_slice_"//licz2//'.stat'
+#endif
+ endif
+ open(istat,file=statname,status="unknown")
+ endif
+
+#ifdef MPI
+ do i=1,scount(me)
+#else
+ do i=1,ntot(islice)
+#endif
+ read(ientout,rec=i,err=101)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+ & eini,efree,rmsdev,(q(k,i),k=1,nQ),iR,ib,iparm
+c write (iout,*) iR,ib,iparm,eini,efree
+ do j=1,2*nres
+ do k=1,3
+ c(k,j)=csingle(k,j)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ iscore=0
+ if (indpdb.gt.0) then
+ call conf_compar(i,.false.,.true.)
+ endif
+ if (bxfile .or.cxfile .or. ensembles.gt.0) write (ientin,rec=i)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+c & potE(i,iparm),-entfac(i),rms_nat,iscore
+ & potE(i,nparmset),-entfac(i),rms_nat,iscore
+c write (iout,'(2i5,3e15.5)') i,me,potE(i,iparm),-entfac(i)
+#ifndef MPI
+ if (cxfile) call cxwrite(ixdrf,csingle,potE(i,nparmset),
+ & -entfac(i),rms_nat,iscore)
+#endif
+ enddo
+ close(ientout,status="delete")
+ close(istat)
+ if (bxfile .or. cxfile .or. ensembles.gt.0) close(ientin)
+#ifdef MPI
+ call MPI_Barrier(WHAM_COMM,IERROR)
+ if (me.ne.Master .or. .not.bxfile .and. .not. cxfile
+ & .and. ensembles.eq.0) return
+ write (iout,*)
+ if (bxfile .or. ensembles.gt.0) then
+ if (nslice.eq.1) then
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//"_par"//licz//".bx"
+ endif
+ else
+ if (.not.separate_parset) then
+ bxname = prefix(:ilen(prefix))//"_slice_"//licz2//".bx"
+ else
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//"par_"//licz//
+ & "_slice_"//licz2//".bx"
+ endif
+ endif
+ open (ientout,file=bxname,status="unknown",
+ & form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Master is creating binary database ",
+ & bxname(:ilen(bxname))
+ endif
+ if (cxfile) then
+ if (nslice.eq.1) then
+ if (.not.separate_parset) then
+ cxname = prefix(:ilen(prefix))//".cx"
+ else
+ cxname = prefix(:ilen(prefix))//"_par"//licz//".cx"
+ endif
+ else
+ if (.not.separate_parset) then
+ cxname = prefix(:ilen(prefix))//
+ & "_slice_"//licz2//".cx"
+ else
+ cxname = prefix(:ilen(prefix))//"_par"//licz//
+ & "_slice_"//licz2//".cx"
+ endif
+ endif
+#if (defined(AIX) && !defined(JUBL))
+ call xdrfopen_(ixdrf,cxname, "w", iret)
+#else
+ call xdrfopen(ixdrf,cxname, "w", iret)
+#endif
+ if (iret.eq.0) then
+ write (iout,*) "Error opening cxfile ",cxname(:ilen(cxname))
+ cxfile=.false.
+ endif
+ endif
+ do j=0,nprocs-1
+ write (liczba,'(bz,i3.3)') j
+ if (separate_parset) then
+ write (licz,'(bz,i3.3)') myparm
+ bxname = prefix(:ilen(prefix))//liczba//"_par"//licz//".bx"
+ else
+ bxname = prefix(:ilen(prefix))//liczba//".bx"
+ endif
+ open (ientin,file=bxname,status="unknown",
+ & form="unformatted",access="direct",recl=lenrec1)
+ write (iout,*) "Master is reading conformations from ",
+ & bxname(:ilen(bxname))
+ iii = 0
+c write (iout,*) "j",j," indstart",indstart(j)," indend",indend(j)
+c call flush(iout)
+ do i=indstart(j),indend(j)
+ iii = iii+1
+ read(ientin,rec=iii,err=101)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+ & eini,efree,rmsdev,iscor
+ if (bxfile .or. ensembles.gt.0) then
+ write (ientout,rec=i)
+ & ((csingle(l,k),l=1,3),k=1,nres),
+ & ((csingle(l,k+nres),l=1,3),k=nnt,nct),
+ & nss,(ihpb(k),jhpb(k),k=1,nss),
+ & eini,efree,rmsdev,iscor
+ endif
+ if(cxfile)call cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
+#ifdef DEBUG
+ do k=1,2*nres
+ do l=1,3
+ c(l,k)=csingle(l,k)
+ enddo
+ enddo
+ call int_from_cart1(.false.)
+ write (iout,'(2i5,3e15.5)') i,iii,eini,efree
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ write (iout,'(16i5)') nss,(ihpb(k),jhpb(k),k=1,nss)
+ write (iout,'(f10.5,i5)') rmsdev,iscor
+#endif
+ enddo ! i
+ write (iout,*) iii," conformations (from",indstart(j)," to",
+ & indend(j),") read from ",
+ & bxname(:ilen(bxname))
+ close (ientin,status="delete")
+ enddo ! j
+ if (bxfile .or. cxfile .or. ensembles.gt.0) close (ientout)
+#if (defined(AIX) && !defined(JUBL))
+ if (cxfile) call xdrfclose_(ixdrf,cxname,iret)
+#else
+ if (cxfile) call xdrfclose(ixdrf,cxname,iret)
+#endif
+#endif
+ return
+ 101 write (iout,*) "Error in scratchfile."
+ call flush(iout)
+ return1
+ end
+c-------------------------------------------------------------------------------
+ subroutine cxwrite(ixdrf,csingle,eini,efree,rmsdev,iscor)
+ implicit none
+ include "DIMENSIONS"
+ include "DIMENSIONS.ZSCOPT"
+ include "DIMENSIONS.FREE"
+ include "DIMENSIONS.COMPAR"
+#ifdef MPI
+ include "mpif.h"
+ integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+ include "COMMON.MPI"
+#endif
+ include "COMMON.CONTROL"
+ include "COMMON.CHAIN"
+ include "COMMON.IOUNITS"
+ include "COMMON.PROTFILES"
+ include "COMMON.NAMES"
+ include "COMMON.VAR"
+ include "COMMON.SBRIDGE"
+ include "COMMON.GEO"
+ include "COMMON.FFIELD"
+ include "COMMON.ENEPS"
+ include "COMMON.LOCAL"
+ include "COMMON.WEIGHTS"
+ include "COMMON.INTERACT"
+ include "COMMON.FREE"
+ include "COMMON.ENERGIES"
+ include "COMMON.COMPAR"
+ include "COMMON.PROT"
+ integer i,j,itmp,iscor,iret,ixdrf
+ double precision rmsdev,efree,eini
+ real*4 csingle(3,maxres2),xoord(3,maxres2+2)
+ real*4 prec
+
+c write (iout,*) "cxwrite"
+c call flush(iout)
+ prec=10000.0
+ do i=1,nres
+ do j=1,3
+ xoord(j,i)=csingle(j,i)
+ enddo
+ enddo
+ do i=nnt,nct
+ do j=1,3
+ xoord(j,nres+i-nnt+1)=csingle(j,i+nres)
+ enddo
+ enddo
+
+ itmp=nres+nct-nnt+1
+
+c write (iout,*) "itmp",itmp
+c call flush(iout)
+#if (defined(AIX) && !defined(JUBL))
+ call xdrf3dfcoord_(ixdrf, xoord, itmp, prec, iret)
+
+c write (iout,*) "xdrf3dfcoord"
+c call flush(iout)
+ call xdrfint_(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint_(ixdrf, ihpb(j), iret)
+ call xdrfint_(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrffloat_(ixdrf,real(eini),iret)
+ call xdrffloat_(ixdrf,real(efree),iret)
+ call xdrffloat_(ixdrf,real(rmsdev),iret)
+ call xdrfint_(ixdrf,iscor,iret)
+#else
+ call xdrf3dfcoord(ixdrf, xoord, itmp, prec, iret)
+
+ call xdrfint(ixdrf, nss, iret)
+ do j=1,nss
+ call xdrfint(ixdrf, ihpb(j), iret)
+ call xdrfint(ixdrf, jhpb(j), iret)
+ enddo
+ call xdrffloat(ixdrf,real(eini),iret)
+ call xdrffloat(ixdrf,real(efree),iret)
+ call xdrffloat(ixdrf,real(rmsdev),iret)
+ call xdrfint(ixdrf,iscor,iret)
+#endif
+
+ return
+ end
+c------------------------------------------------------------------------------
+ logical function conf_check(ii,iprint)
+ implicit none
+ include "DIMENSIONS"
+ include "DIMENSIONS.ZSCOPT"
+ include "DIMENSIONS.FREE"
+#ifdef MPI
+ include "mpif.h"
+ integer IERROR,ERRCODE,STATUS(MPI_STATUS_SIZE)
+ include "COMMON.MPI"
+#endif
+ include "COMMON.CHAIN"
+ include "COMMON.IOUNITS"
+ include "COMMON.PROTFILES"
+ include "COMMON.NAMES"
+ include "COMMON.VAR"
+ include "COMMON.SBRIDGE"
+ include "COMMON.GEO"
+ include "COMMON.FFIELD"
+ include "COMMON.ENEPS"
+ include "COMMON.LOCAL"
+ include "COMMON.WEIGHTS"
+ include "COMMON.INTERACT"
+ include "COMMON.FREE"
+ include "COMMON.ENERGIES"
+ include "COMMON.CONTROL"
+ include "COMMON.TORCNSTR"
+ integer j,k,l,ii,itj,iprint
+ if (.not.check_conf) then
+ conf_check=.true.
+ return
+ endif
+ call int_from_cart1(.false.)
+ do j=nnt+1,nct
+ if (vbld(j).lt.2.0d0 .or. vbld(j).gt.5.0d0) then
+ if (iprint.gt.0)
+ & write (iout,*) "Bad CA-CA bond length",j," ",vbld(j),
+ & " for conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*)
+ & "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ enddo
+ do j=nnt,nct
+ itj=itype(j)
+ if (itype(j).ne.10 .and. (vbld(nres+j)-dsc(itj)).gt.2.0d0) then
+ if (iprint.gt.0)
+ & write (iout,*) "Bad CA-SC bond length",j," ",vbld(nres+j),
+ & " for conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*)
+ & "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ enddo
+ do j=3,nres
+ if (theta(j).le.0.0d0) then
+ if (iprint.gt.0)
+ & write (iout,*) "Zero theta angle(s) in conformation",ii
+ if (iprint.gt.1) then
+ write (iout,*) "The Cartesian geometry is:"
+ write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
+ write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
+ write (iout,*) "The internal geometry is:"
+ write (iout,'(8f10.4)') (vbld(k),k=nnt+1,nct)
+ write (iout,'(8f10.4)') (vbld(k),k=nres+nnt,nres+nct)
+ write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
+ write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
+ write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
+ write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
+ endif
+ if (iprint.gt.0) write (iout,*)
+ & "This conformation WILL NOT be added to the database."
+ conf_check=.false.
+ return
+ endif
+ if (theta(j).gt.179.97*deg2rad) theta(j)=179.97*deg2rad
+ enddo
+ conf_check=.true.
+c write (iout,*) "conf_check passed",ii
+ return
+ end
projects / unres.git / blobdiff