--- /dev/null
+cc---------------------------------
+ subroutine refresh_bank(ntrial)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'mpif.h'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CHAIN'
+ include 'COMMON.VAR'
+ include 'COMMON.CONTROL'
+ character chacc
+ integer iaccn
+ double precision l_diff(mxio),denep
+
+ do i=0,mxmv
+ do j=1,3
+ nstatnx_tot(i,j)=nstatnx_tot(i,j)+nstatnx(i,j)
+ nstatnx(i,j)=0
+ enddo
+ enddo
+
+c loop over all newly obtained conformations
+ do n=1,ntrial
+ chacc=' '
+ iaccn=0
+ nstatnx(movernx(n),1)=nstatnx(movernx(n),1)+1
+cccccccccccccccccccccccccccccccccccccccccccc
+cjlee
+ if(iref.ne.0) then
+ if(rmsn(n).gt.rmscut.or.pncn(n).lt.pnccut) goto 100
+ endif
+cjlee
+ if(etot(n).gt.ebmax) goto 100
+c Find the conformation closest to the conformation n in the bank
+ difmin=9.d9
+ do m=1,nbank
+ call get_diff12(dihang(1,1,1,n),bvar(1,1,1,m),l_diff(m))
+ if(l_diff(m).lt.difmin) then
+ difmin=l_diff(m)
+ idmin=m
+ endif
+ enddo
+
+ if(difmin.lt.cutdif) then
+c n is redundant to idmin
+ if(etot(n).lt.bene(idmin)) then
+ if(etot(n).lt.bene(idmin)-0.01d0) then
+ ibank(idmin)=0
+ jbank(idmin)=0
+ endif
+ denep=bene(idmin)-etot(n)
+ call replace_bvar(idmin,n)
+crc Update dij
+ do i1=1,nbank
+ if (i1.ne.idmin) then
+ dij(i1,idmin)=l_diff(i1)
+ dij(idmin,i1)=l_diff(i1)
+ endif
+ enddo
+ chacc='c'
+ iaccn=idmin
+ nstatnx(movernx(n),2)=nstatnx(movernx(n),2)+1
+ if(idmin.eq.ibmax) call find_max
+ endif
+ else
+c got new conformation
+ del_ene=0.0d0
+ if(ebmax-ebmin.gt.del_ene) then
+ denep=ebmax-etot(n)
+ call replace_bvar(ibmax,n)
+crc Update dij
+ do i1=1,nbank
+ if (i1.ne.ibmax) then
+ dij(i1,ibmax)=l_diff(i1)
+ dij(ibmax,i1)=l_diff(i1)
+ endif
+ enddo
+ chacc='f'
+ iaccn=ibmax
+ nstatnx(movernx(n),3)=nstatnx(movernx(n),3)+1
+ ibank(ibmax)=0
+ jbank(ibmax)=0
+ call find_max
+ else
+ if(del_ene.lt.0.0001) then
+ write (iout,*) 'ERROR in refresh_bank: '
+ write (iout,*) 'ebmax: ',ebmax
+ write (iout,*) 'ebmin: ',ebmin
+ write (iout,*) 'del_ene: ',del_ene
+crc call mpi_abort(mpi_comm_world,ierror,ierrcode)
+ endif
+cjp nbmax is never defined so condition below is always false
+c if(nbank.lt.nbmax) then
+c nbank=nbank+1
+c call replace_bvar(nbank,n)
+c ibank(nbank)=0
+c jbank(nbank)=0
+c else
+ call replace_bvar(ibmax,n)
+ ibank(ibmax)=0
+ jbank(ibmax)=0
+ call find_max
+c endif
+ endif
+ endif
+cccccccccccccccccccccccccccccccccccccccccccc
+ 100 continue
+ if (iaccn.eq.0) then
+ if (iref.eq.0) then
+ write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a4,i3,i4,3i5)')
+ & indb(n,2),' e ',indb(n,3),indb(n,1),' etot ',etot(n),' mv ',
+ & indb(n,5),indb(n,4),indb(n,7),indb(n,8),indb(n,9)
+ else
+ write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a4,i3,i4,3i5
+ & ,a5,0pf4.1,a5,f3.0)')
+ & indb(n,2),' e ',indb(n,3),indb(n,1),' etot ',etot(n),' mv ',
+ & indb(n,5),indb(n,4),indb(n,7),indb(n,8),indb(n,9),
+ & ' rms ',rmsn(n),' %NC ',pncn(n)*100
+ endif
+ else
+ if (iref.eq.0) then
+ write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a4,i3,i4,3i5,
+ & 1x,a1,i4,0pf8.1,0pf8.1)')
+ & indb(n,2),' e ',indb(n,3),indb(n,1),' etot ',etot(n),' mv ',
+ & indb(n,5),indb(n,4),indb(n,7),indb(n,8),indb(n,9),
+ & chacc,iaccn,difmin,denep
+ else
+ write(iout,'(i3,a3,i4,i5,a6,1pe12.4,a4,i3,i4,3i5,a5,
+ & 0pf4.1,a5,f3.0,1x,a1,i4,0pf8.1,0pf8.1)')
+ & indb(n,2),' e ',indb(n,3),indb(n,1),' etot ',etot(n),' mv ',
+ & indb(n,5),indb(n,4),indb(n,7),indb(n,8),indb(n,9),
+ & ' rms ',rmsn(n),' %NC ',pncn(n)*100,
+ & chacc,iaccn,difmin,denep
+ endif
+ endif
+ enddo
+c end of loop over all newly obtained conformations
+ do i=0,mxmv
+ if(nstatnx(i,1).ne.0) then
+ if (i.le.9) then
+ write(iout,'(a4,i1,a7,i4,a7,i4,a5,i4,a5,f5.1)')
+ & '## N',i,' total=',nstatnx(i,1),
+ & ' close=',nstatnx(i,2),' far=',nstatnx(i,3),
+ & ' %acc',(nstatnx(i,2)+nstatnx(i,3))*100.0/nstatnx(i,1)
+ else
+ write(iout,'(a3,i2,a7,i4,a7,i4,a5,i4,a5,f5.1)')
+ & '##N',i,' total=',nstatnx(i,1),
+ & ' close=',nstatnx(i,2),' far=',nstatnx(i,3),
+ & ' %acc',(nstatnx(i,2)+nstatnx(i,3))*100.0/nstatnx(i,1)
+ endif
+ else
+ if (i.le.9) then
+ write(iout,'(a4,i1,a7,i4,a7,i4,a5,i4,a5,f5.1)')
+ & '## N',i,' total=',nstatnx(i,1),
+ & ' close=',nstatnx(i,2),' far=',nstatnx(i,3),
+ & ' %acc',0.0
+ else
+ write(iout,'(a3,i2,a7,i4,a7,i4,a5,i4,a5,f5.1)')
+ & '##N',i,' total=',nstatnx(i,1),
+ & ' close=',nstatnx(i,2),' far=',nstatnx(i,3),
+ & ' %acc',0.0
+ endif
+ endif
+ enddo
+ call flush(iout)
+crc Update dij
+crc moved up, saves some get_diff12 calls
+crc
+crc do i1=1,nbank-1
+crc do i2=i1+1,nbank
+crc if(jbank(i1).eq.0.or.jbank(i2).eq.0) then
+crc call get_diff12(bvar(1,1,1,i1),bvar(1,1,1,i2),diff)
+crc dij(i1,i2)=diff
+crc dij(i2,i1)=diff
+crc endif
+crc enddo
+crc enddo
+
+ do i=1,nbank
+ jbank(i)=1
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine replace_bvar(iold,inew)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'mpif.h'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.CONTROL'
+ include 'COMMON.SBRIDGE'
+
+ if (iold.gt.mxio .or. iold.lt.1 .or. inew.gt.mxio .or. inew.lt.1)
+ & then
+ write (iout,*) 'Dimension ERROR in REPLACE_BVAR: IOLD',iold,
+ & ' INEW',inew
+ call mpi_abort(mpi_comm_world,ierror,ierrcode)
+ endif
+ do k=1,numch
+ do j=2,nres-1
+ do i=1,4
+ bvar(i,j,k,iold)=dihang(i,j,k,inew)
+ enddo
+ enddo
+ enddo
+ bene(iold)=etot(inew)
+ brmsn(iold)=rmsn(inew)
+ bpncn(iold)=pncn(inew)
+
+ if(bene(iold).lt.ebmin) then
+ ebmin=bene(iold)
+ ibmin=iold
+ endif
+
+ if(vdisulf) then
+ bvar_nss(iold)=nss_out(inew)
+cd write(iout,*) 'SS BANK',iold,bvar_nss(iold)
+ do i=1,bvar_nss(iold)
+ bvar_ss(1,i,iold)=iss_out(i,inew)
+ bvar_ss(2,i,iold)=jss_out(i,inew)
+cd write(iout,*) 'SS',bvar_ss(1,i,iold)-nres,
+cd & bvar_ss(2,i,iold)-nres
+ enddo
+
+ bvar_ns(iold)=ns-2*bvar_nss(iold)
+cd write(iout,*) 'CYS #free ', bvar_ns(iold)
+ k=0
+ do i=1,ns
+ j=1
+ do while( iss(i).ne.iss_out(j,inew)-nres .and.
+ & iss(i).ne.jss_out(j,inew)-nres .and.
+ & j.le.nss_out(inew))
+ j=j+1
+ enddo
+ if (j.gt.nss_out(inew)) then
+ k=k+1
+ bvar_s(k,iold)=iss(i)
+ endif
+ enddo
+cd write(iout,*) 'CYS free',(bvar_s(k,iold),k=1,bvar_ns(iold))
+ endif
+
+ return
+ end
+c---------------------------------------
+ subroutine write_rbank(jlee,adif,nft)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+
+ open(icsa_rbank,file=csa_rbank,status="unknown")
+ write (icsa_rbank,900) jlee,nbank,nstep,nft,icycle,adif
+ do k=1,nbank
+ write (icsa_rbank,952) k,rene(k),rrmsn(k),rpncn(k)
+ do j=1,numch
+ do l=2,nres-1
+ write (icsa_rbank,850) (rad2deg*rvar(i,l,j,k),i=1,4)
+ enddo
+ enddo
+ enddo
+ close(icsa_rbank)
+
+ 850 format (10f8.3)
+ 900 format (1x,"jlee =",i3,3x,"nbank =",i4,3x,"nstep =",
+ & i8,i10,i2,f15.5)
+ 952 format (1x,'#',i4,' total E ',1pe14.5,' rmsd from N ',0pf8.3
+ & ,' %NC ',0pf5.2)
+
+ return
+ end
+c---------------------------------------
+ subroutine read_rbank(jlee,adif)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'mpif.h'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ include 'COMMON.SETUP'
+ character*80 karta
+
+ open(icsa_rbank,file=csa_rbank,status="old")
+ read (icsa_rbank,901) jleer,nbankr,nstepr,nftr,icycler,adif
+ print *,jleer,nbankr,nstepr,nftr,icycler,adif
+c print *, 'adif from read_rbank ',adif
+ if(nbankr.ne.nbank) then
+ write (iout,*) 'ERROR in READ_BANK: NBANKR',nbankr,
+ & ' NBANK',nbank
+ call mpi_abort(mpi_comm_world,ierror,ierrcode)
+ endif
+ if(jleer.ne.jlee) then
+ write (iout,*) 'ERROR in READ_BANK: JLEER',jleer,
+ & ' JLEE',jlee
+ call mpi_abort(mpi_comm_world,ierror,ierrcode)
+ endif
+
+ kk=0
+ do k=1,nbankr
+ read (icsa_rbank,'(a80)') karta
+ write(iout,*) "READ_RBANK: kk=",kk
+ write(iout,*) karta
+c if (index(karta,"*").gt.0) then
+c write (iout,*) "***** Stars in bankr ***** k=",k,
+c & " skipped"
+c do j=1,numch
+c do l=2,nres-1
+c read (30,850) (rdummy,i=1,4)
+c enddo
+c enddo
+c else
+ kk=kk+1
+ call reada(karta,"total E",rene(kk),1.0d20)
+ call reada(karta,"rmsd from N",rrmsn(kk),0.0d0)
+ call reada(karta,"%NC",rpncn(kk),0.0d0)
+ write(iout,*)"total E",bene(kk),"rmsd from N",brmsn(kk),
+ & "%NC",bpncn(kk),ibank(kk)
+c read (icsa_rbank,953) kdummy,rene(kk),rrmsn(kk),rpncn(kk)
+ do j=1,numch
+ do l=2,nres-1
+ read (icsa_rbank,850) (rvar(i,l,j,kk),i=1,4)
+c write (iout,850) (rvar(i,l,j,kk),i=1,4)
+ do i=1,4
+ rvar(i,l,j,kk)=deg2rad*rvar(i,l,j,kk)
+ enddo
+ enddo
+ enddo
+c endif
+ enddo
+cd write (*,*) "read_rbank ******************* kk",kk,
+cd & "nbankr",nbankr
+ if (kk.lt.nbankr) nbankr=kk
+cd do kk=1,nbankr
+cd print *,"kk=",kk
+cd do j=1,numch
+cd do l=2,nres-1
+cd write (*,850) (rvar(i,l,j,kk),i=1,4)
+cd enddo
+cd enddo
+cd enddo
+ close(icsa_rbank)
+
+ 850 format (10f8.3)
+ 901 format (1x,6x,i3,3x,7x,i4,3x,7x,i8,i10,i2,f15.5)
+ 953 format (1x,1x,i4,9x,f12.3,13x,f8.3,5x,f5.2)
+
+ return
+ end
+c---------------------------------------
+ subroutine write_bank(jlee,nft)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.CONTROL'
+ character*7 chtmp
+ character*40 chfrm
+ external ilen
+
+ open(icsa_bank,file=csa_bank,status="unknown")
+ write (icsa_bank,900) jlee,nbank,nstep,nft,icycle,cutdif
+ write (icsa_bank,902) nglob_csa, eglob_csa
+ open (igeom,file=intname,status='UNKNOWN')
+ do k=1,nbank
+ write (icsa_bank,952) k,bene(k),brmsn(k),bpncn(k),ibank(k)
+ if (vdisulf) write (icsa_bank,'(101i4)')
+ & bvar_nss(k),((bvar_ss(j,i,k),j=1,2),i=1,bvar_nss(k))
+ do j=1,numch
+ do l=2,nres-1
+ write (icsa_bank,850) (rad2deg*bvar(i,l,j,k),i=1,4)
+ enddo
+ enddo
+ if (bvar_nss(k).le.9) then
+ write (igeom,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & bvar_nss(k),(bvar_ss(1,i,k),bvar_ss(2,i,k),i=1,bvar_nss(k))
+ else
+ write (igeom,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & bvar_nss(k),(bvar_ss(1,i,k),bvar_ss(2,i,k),i=1,9)
+ write (igeom,'(3X,11(1X,2I3))') (bvar_ss(1,i,k),
+ & bvar_ss(2,i,k),i=10,bvar_nss(k))
+ endif
+ write (igeom,200) (rad2deg*bvar(1,i,1,k),i=2,nres-1)
+ write (igeom,200) (rad2deg*bvar(2,i,1,k),i=2,nres-2)
+ write (igeom,200) (rad2deg*bvar(3,i,1,k),i=2,nres-1)
+ write (igeom,200) (rad2deg*bvar(4,i,1,k),i=2,nres-1)
+ enddo
+ close(icsa_bank)
+ close(igeom)
+
+ if (nstep/200.gt.ilastnstep) then
+
+ ilastnstep=(ilastnstep+1)*1.5
+ write(chfrm,'(a2,i1,a1)') '(i',int(dlog10(dble(nstep))+1),')'
+ write(chtmp,chfrm) nstep
+ open(icsa_int,file=prefix(:ilen(prefix))
+ & //'_'//chtmp(:ilen(chtmp))//'.int',status='UNKNOWN')
+ do k=1,nbank
+ if (bvar_nss(k).le.9) then
+ write (icsa_int,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & bvar_nss(k),(bvar_ss(1,i,k),bvar_ss(2,i,k),i=1,bvar_nss(k))
+ else
+ write (icsa_int,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & bvar_nss(k),(bvar_ss(1,i,k),bvar_ss(2,i,k),i=1,9)
+ write (icsa_int,'(3X,11(1X,2I3))') (bvar_ss(1,i,k),
+ & bvar_ss(2,i,k),i=10,bvar_nss(k))
+ endif
+ write (icsa_int,200) (rad2deg*bvar(1,i,1,k),i=2,nres-1)
+ write (icsa_int,200) (rad2deg*bvar(2,i,1,k),i=2,nres-2)
+ write (icsa_int,200) (rad2deg*bvar(3,i,1,k),i=2,nres-1)
+ write (icsa_int,200) (rad2deg*bvar(4,i,1,k),i=2,nres-1)
+ enddo
+ close(icsa_int)
+ endif
+
+
+ 200 format (8f10.4)
+ 850 format (10f8.3)
+ 900 format (1x,"jlee =",i3,3x,"nbank =",i4,3x,"nstep =",
+ & i8,i10,i2,f15.5)
+ 902 format (1x,'nglob_csa =',i4,' eglob_csa =',1pe14.5)
+ 952 format (1x,'#',i4,' total E ',1pe14.5,' rmsd from N ',0pf8.3,
+ & ' %NC ',0pf5.2,i5)
+
+ return
+ end
+c---------------------------------------
+ subroutine write_bank_reminimized(jlee,nft)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ include 'COMMON.SBRIDGE'
+
+ open(icsa_bank_reminimized,file=csa_bank_reminimized,
+ & status="unknown")
+ write (icsa_bank_reminimized,900)
+ & jlee,nbank,nstep,nft,icycle,cutdif
+ open (igeom,file=intname,status='UNKNOWN')
+ do k=1,nbank
+ write (icsa_bank_reminimized,952) k,bene(k),brmsn(k),
+ & bpncn(k),ibank(k)
+ do j=1,numch
+ do l=2,nres-1
+ write (icsa_bank_reminimized,850) (rad2deg*bvar(i,l,j,k),i=1,4)
+ enddo
+ enddo
+ if (nss.le.9) then
+ write (igeom,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & nss,(ihpb(i),jhpb(i),i=1,nss)
+ else
+ write (igeom,'(I5,F10.3,I2,9(1X,2I3))') k,bene(k),
+ & nss,(ihpb(i),jhpb(i),i=1,9)
+ write (igeom,'(3X,11(1X,2I3))') (ihpb(i),jhpb(i),i=10,nss)
+ endif
+ write (igeom,200) (rad2deg*bvar(1,i,1,k),i=2,nres-1)
+ write (igeom,200) (rad2deg*bvar(2,i,1,k),i=2,nres-2)
+ write (igeom,200) (rad2deg*bvar(3,i,1,k),i=2,nres-1)
+ write (igeom,200) (rad2deg*bvar(4,i,1,k),i=2,nres-1)
+ enddo
+ close(icsa_bank_reminimized)
+ close(igeom)
+
+ 200 format (8f10.4)
+ 850 format (10f8.3)
+ 900 format (1x,"jlee =",i3,3x,"nbank =",i4,3x,"nstep =",
+ & i8,i10,i2,f15.5)
+ 952 format (1x,'#',i4,' total E ',1pe14.5,' rmsd from N ',0pf8.3
+ & ,' %NC ',0pf5.2,i5)
+
+ return
+ end
+c---------------------------------
+ subroutine read_bank(jlee,nft,cutdifr)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+ include 'COMMON.CONTROL'
+ include 'COMMON.SBRIDGE'
+ character*80 karta
+ integer ilen
+ external ilen
+
+ open(icsa_bank,file=csa_bank,status="old")
+ read (icsa_bank,901) jlee,nbank,nstep,nft,icycle,cutdifr
+ read (icsa_bank,902) nglob_csa, eglob_csa
+c if(jleer.ne.jlee) then
+c write (iout,*) 'ERROR in READ_BANK: JLEER',jleer,
+c & ' JLEE',jlee
+c call mpi_abort(mpi_comm_world,ierror,ierrcode)
+c endif
+
+ kk=0
+ do k=1,nbank
+ read (icsa_bank,'(a80)') karta
+ write(iout,*) "READ_BANK: kk=",kk
+ write(iout,*) karta
+c if (index(karta,"*").gt.0) then
+c write (iout,*) "***** Stars in bank ***** k=",k,
+c & " skipped"
+c do j=1,numch
+c do l=2,nres-1
+c read (33,850) (rdummy,i=1,4)
+c enddo
+c enddo
+c else
+ kk=kk+1
+ call reada(karta,"total E",bene(kk),1.0d20)
+ call reada(karta,"rmsd from N",brmsn(kk),0.0d0)
+ call reada(karta,"%NC",bpncn(kk),0.0d0)
+ read (karta(ilen(karta)-1:),*,end=111,err=111) ibank(kk)
+ goto 112
+ 111 ibank(kk)=0
+ 112 continue
+ write(iout,*)"total E",bene(kk),"rmsd from N",brmsn(kk),
+ & "%NC",bpncn(kk),ibank(kk)
+c read (icsa_bank,953) kdummy,bene(k),brmsn(k),bpncn(k),ibank(k)
+ if (vdisulf) then
+ read (icsa_bank,'(101i4)')
+ & bvar_nss(kk),((bvar_ss(j,i,kk),j=1,2),i=1,bvar_nss(kk))
+ bvar_ns(kk)=ns-2*bvar_nss(kk)
+ write(iout,*) 'read SSBOND',bvar_nss(kk),
+ & ((bvar_ss(j,i,kk),j=1,2),i=1,bvar_nss(kk))
+cd write(iout,*) 'read CYS #free ', bvar_ns(kk)
+ l=0
+ do i=1,ns
+ j=1
+ do while( iss(i).ne.bvar_ss(1,j,kk)-nres .and.
+ & iss(i).ne.bvar_ss(2,j,kk)-nres .and.
+ & j.le.bvar_nss(kk))
+ j=j+1
+ enddo
+ if (j.gt.bvar_nss(kk)) then
+ l=l+1
+ bvar_s(l,kk)=iss(i)
+ endif
+ enddo
+cd write(iout,*)'read CYS free',(bvar_s(l,kk),l=1,bvar_ns(kk))
+ endif
+ do j=1,numch
+ do l=2,nres-1
+ read (icsa_bank,850) (bvar(i,l,j,kk),i=1,4)
+c write (iout,850) (bvar(i,l,j,kk),i=1,4)
+ do i=1,4
+ bvar(i,l,j,kk)=deg2rad*bvar(i,l,j,kk)
+ enddo ! l
+ enddo ! l
+ enddo ! j
+c endif
+ enddo ! k
+
+ if (kk.lt.nbank) nbank=kk
+cd write (*,*) "read_bank ******************* kk",kk,
+cd & "nbank",nbank
+cd do kk=1,nbank
+cd print *,"kk=",kk
+cd do j=1,numch
+cd do l=2,nres-1
+cd write (*,850) (bvar(i,l,j,kk),i=1,4)
+cd enddo
+cd enddo
+cd enddo
+
+c do k=1,nbank
+c read (33,953) kdummy,bene(k),brmsn(k),bpncn(k),ibank(k)
+c do j=1,numch
+c do l=2,nres-1
+c read (33,850) (bvar(i,l,j,k),i=1,4)
+c do i=1,4
+c bvar(i,l,j,k)=deg2rad*bvar(i,l,j,k)
+c enddo
+c enddo
+c enddo
+c enddo
+ close(icsa_bank)
+
+ 850 format (10f8.3)
+ 952 format (1x,'#',i4,' total E ',f12.3,' rmsd from N ',f8.3,i5)
+ 901 format (1x,6x,i3,3x,7x,i4,3x,7x,i8,i10,i2,f15.5)
+ 902 format (1x,11x,i4,12x,1pe14.5)
+ 953 format (1x,1x,i4,9x,f12.3,13x,f8.3,5x,f5.2,i5)
+
+ return
+ end
+c---------------------------------------
+ subroutine write_bank1(jlee)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+ include 'COMMON.GEO'
+
+#if defined(AIX) || defined(PGI)
+ open(icsa_bank1,file=csa_bank1,position="append")
+#else
+ open(icsa_bank1,file=csa_bank1,access="append")
+#endif
+ write (icsa_bank1,900) jlee,nbank,nstep,cutdif
+ do k=1,nbank
+ write (icsa_bank1,952) k,bene(k),brmsn(k),bpncn(k),ibank(k)
+ do j=1,numch
+ do l=2,nres-1
+ write (icsa_bank1,850) (rad2deg*bvar(i,l,j,k),i=1,4)
+ enddo
+ enddo
+ enddo
+ close(icsa_bank1)
+ 850 format (10f8.3)
+ 900 format (4x,"jlee =",i5,3x,"nbank =",i5,3x,"nstep =",i10,f15.5)
+ 952 format (1x,'#',i4,' total E ',1pe14.5,' rmsd from N ',0pf8.3
+ & ,' %NC ',0pf5.2,i5)
+
+ return
+ end
+c---------------------------------
+ subroutine save_is(ind)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'mpif.h'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.CHAIN'
+
+ index=nbank+ind
+c print *, "nbank,ind,index,is(ind) ",nbank,ind,index,is(ind)
+ if (index.gt.mxio .or. index.lt.1 .or.
+ & is(ind).gt.mxio .or. is(ind).lt.1) then
+ write (iout,*) 'Dimension ERROR in SAVE_IS: INDEX',index,
+ & ' IND',ind,' IS',is(ind)
+ call mpi_abort(mpi_comm_world,ierror,ierrcode)
+ endif
+ do k=1,numch
+ do j=2,nres-1
+ do i=1,4
+ bvar(i,j,k,index)=bvar(i,j,k,is(ind))
+ enddo
+ enddo
+ enddo
+ bene(index)=bene(is(ind))
+ ibank(is(ind))=1
+
+ return
+ end
+c---------------------------------
+ subroutine select_is(n,ifar,idum)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ dimension itag(mxio),adiff(mxio)
+
+ iuse=0
+ do i=1,nbank
+ if(ibank(i).eq.0) then
+ iuse=iuse+1
+ itag(iuse)=i
+ endif
+ enddo
+ iusesv=iuse
+
+ if(iuse.eq.0) then
+ icycle=icycle+1
+ do i=1,nbank
+ if(ibank(i).eq.2) then
+ ibank(i)=1
+ else
+ ibank(i)=0
+ endif
+ enddo
+ imade=0
+ call get_is(idum,ifar,n,imade,0)
+ctest3 call get_is_max(idum,ifar,n,imade,0)
+ else if(iuse.eq.n) then
+ do i=1,iuse
+ is(i)=itag(i)
+ call save_is(i)
+ enddo
+ else if(iuse.lt.n) then
+c if(icycle.eq.0) then
+c do i=1,n
+c ind=mod(i-1,iuse)+1
+c is(i)=itag(ind)
+c call save_is(i)
+c enddo
+c else
+c endif
+ do i=1,iuse
+ is(i)=itag(i)
+ call save_is(i)
+ enddo
+ imade=iuse
+c call get_is_ran(idum,n,imade,1)
+ call get_is(idum,ifar,n,imade,1)
+ctest3 call get_is_max(idum,ifar,n,imade,1)
+c if(iusesv.le.n/10) then
+ if(iusesv.le.0) then
+ icycle=icycle+1
+ do i=1,nbank
+c if(ibank(i).eq.2) then
+c ibank(i)=1
+ if(ibank(i).ge.2) then
+ ibank(i)=ibank(i)-1
+ else
+ ibank(i)=0
+ endif
+ enddo
+ endif
+ else
+ imade=0
+ call get_is(idum,ifar,n,imade,0)
+ctest3 call get_is_max(idum,ifar,n,imade,0)
+ endif
+ iuse=iusesv
+
+ return
+ end
+c---------------------------------
+ subroutine get_is_ran(idum,n,imade,k)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ real ran1,ran2
+ dimension itag(mxio),adiff(mxio)
+
+ do j=imade+1,n
+ iuse=0
+ do i=1,nbank
+ if(ibank(i).eq.k) then
+ iuse=iuse+1
+ itag(iuse)=i
+ endif
+ enddo
+ iran=iuse* ran1(idum)+1
+ is(j)=itag(iran)
+ call save_is(j)
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine get_is(idum,ifar,n,imade,k)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ real ran1,ran2
+ dimension itag(mxio),adiff(mxio)
+
+ iuse=0
+ do i=1,nbank
+ if(ibank(i).eq.k) then
+ iuse=iuse+1
+ itag(iuse)=i
+ endif
+ enddo
+ iran=iuse* ran1(idum)+1
+ imade=imade+1
+ is(imade)=itag(iran)
+ call save_is(imade)
+
+ do i=imade+1,ifar-1
+ if(icycle.eq.-1) then
+ call select_iseed_max(i,k)
+ else
+ call select_iseed_min(i,k)
+ctest4 call select_iseed_max(i,k)
+ endif
+ call save_is(i)
+ enddo
+
+ do i=ifar,n
+ call select_iseed_far(i,k)
+ call save_is(i)
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine select_iseed_max(imade1,ik)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ dimension itag(mxio),adiff(mxio)
+
+ iuse=0
+ avedif=0.d0
+ difmax=0.d0
+ do n=1,nbank
+ if(ibank(n).eq.ik) then
+ iuse=iuse+1
+ diffmn=9.d190
+ do imade=1,imade1-1
+c m=nbank+imade
+c call get_diff12(bvar(1,1,1,n),bvar(1,1,1,m),diff,idiff)
+ m=is(imade)
+ diff=dij(n,m)
+ if(diff.lt.diffmn) diffmn=diff
+ enddo
+ if(diffmn.gt.difmax) difmax=diffmn
+ adiff(iuse)=diffmn
+ itag(iuse)=n
+ avedif=avedif+diffmn
+ endif
+ enddo
+
+ avedif=avedif/iuse
+c avedif=(avedif+difmax)/2
+ emax=-9.d190
+ do i=1,iuse
+ if(adiff(i).ge.avedif) then
+ itagi=itag(i)
+ benei=bene(itagi)
+ if(benei.gt.emax) then
+ emax=benei
+ is(imade1)=itagi
+ endif
+ endif
+ enddo
+
+ if(ik.eq.0) iuse=iuse-1
+
+ return
+ end
+c---------------------------------
+ subroutine select_iseed_min(imade1,ik)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ dimension itag(mxio),adiff(mxio)
+
+ iuse=0
+ avedif=0.d0
+ difmax=0.d0
+ do n=1,nbank
+ if(ibank(n).eq.ik) then
+ iuse=iuse+1
+ diffmn=9.d190
+ do imade=1,imade1-1
+c m=nbank+imade
+c call get_diff12(bvar(1,1,1,n),bvar(1,1,1,m),diff,idiff)
+ m=is(imade)
+ diff=dij(n,m)
+ if(diff.lt.diffmn) diffmn=diff
+ enddo
+ if(diffmn.gt.difmax) difmax=diffmn
+ adiff(iuse)=diffmn
+ itag(iuse)=n
+ avedif=avedif+diffmn
+ endif
+ enddo
+
+ avedif=avedif/iuse
+c avedif=(avedif+difmax)/2
+ emin=9.d190
+ do i=1,iuse
+c print *,"i, adiff(i),avedif : ",i,adiff(i),avedif
+ if(adiff(i).ge.avedif) then
+ itagi=itag(i)
+ benei=bene(itagi)
+c print *,"i, benei,emin : ",i,benei,emin
+ if(benei.lt.emin) then
+ emin=benei
+ is(imade1)=itagi
+ endif
+ endif
+ enddo
+
+ if(ik.eq.0) iuse=iuse-1
+
+c print *, "exiting select_iseed_min",is(imade1)
+
+ return
+ end
+c---------------------------------
+ subroutine select_iseed_far(imade1,ik)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+
+ dmax=-9.d190
+ do n=1,nbank
+ if(ibank(n).eq.ik) then
+ diffmn=9.d190
+ do imade=1,imade1-1
+c m=nbank+imade
+c call get_diff12(bvar(1,1,1,n),bvar(1,1,1,m),diff,idiff)
+ m=is(imade)
+ diff=dij(n,m)
+ if(diff.lt.diffmn) diffmn=diff
+ enddo
+ endif
+ if(diffmn.gt.dmax) then
+ dmax=diffmn
+ is(imade1)=n
+ endif
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine find_min
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+
+ ebmin=9.d190
+
+ do i=1,nbank
+ benei=bene(i)
+ if(benei.lt.ebmin) then
+ ebmin=benei
+ ibmin=i
+ endif
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine write_csa_pdb(var,ene,nft,ik,iw_pdb)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ include 'COMMON.VAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.MINIM'
+ include 'COMMON.SETUP'
+ include 'COMMON.GEO'
+ include 'COMMON.CHAIN'
+ include 'COMMON.LOCAL'
+ include 'COMMON.INTERACT'
+ include 'COMMON.NAMES'
+ include 'COMMON.SBRIDGE'
+ integer lenpre,lenpot,ilen
+ external ilen
+ dimension var(maxvar)
+ character*50 titelloc
+ character*3 zahl
+
+ nmin_csa=nmin_csa+1
+ if(ene.lt.eglob_csa) then
+ eglob_csa=ene
+ nglob_csa=nglob_csa+1
+ call numstr(nglob_csa,zahl)
+
+ call var_to_geom(nvar,var)
+ call chainbuild
+ call secondary2(.false.)
+
+ lenpre=ilen(prefix)
+ open(icsa_pdb,file=prefix(:lenpre)//'@'//zahl//'.pdb')
+
+ if (iw_pdb.eq.1) then
+ write(titelloc,'(a2,i3,a3,i9,a3,i6)')
+ & 'GM',nglob_csa,' e ',nft,' m ',nmin_csa
+ else
+ write(titelloc,'(a2,i3,a3,i9,a3,i6,a5,f5.2,a5,f5.1)')
+ & 'GM',nglob_csa,' e ',nft,' m ',nmin_csa,' rms '
+ & ,rmsn(ik),' %NC ',pncn(ik)*100
+ endif
+ call pdbout(eglob_csa,titelloc,icsa_pdb)
+ close(icsa_pdb)
+ endif
+
+ return
+ end
+c---------------------------------
+ subroutine find_max
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+
+ ebmax=-9.d190
+
+ do i=1,nbank
+ benei=bene(i)
+ if(benei.gt.ebmax) then
+ ebmax=benei
+ ibmax=i
+ endif
+ enddo
+
+ return
+ end
+c---------------------------------
+ subroutine get_diff
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+
+ tdiff=0.d0
+ difmin=9.d190
+ do i1=1,nbank-1
+ do i2=i1+1,nbank
+ if(jbank(i1).eq.0.or.jbank(i2).eq.0) then
+ call get_diff12(bvar(1,1,1,i1),bvar(1,1,1,i2),diff)
+ dij(i1,i2)=diff
+ dij(i2,i1)=diff
+ else
+ diff=dij(i1,i2)
+ endif
+ tdiff=tdiff+diff
+ if(diff.lt.difmin) difmin=diff
+ enddo
+ dij(i1,i1)=0.0
+ enddo
+
+ do i=1,nbank
+ jbank(i)=1
+ enddo
+
+ avedif=tdiff/nbank/(nbank-1)*2
+
+ return
+ end
+c---------------------------------
+ subroutine estimate_cutdif(adif,xct,cutdifr)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+
+ ctdif1=adif/cut2
+
+ exponent = cutdifr*cut1/adif
+ exponent = dlog(exponent)/dlog(xct)
+
+ nexp=exponent+0.25
+ cutdif= adif/cut1*xct**nexp
+ if(cutdif.lt.ctdif1) cutdif=ctdif1
+
+ return
+ end
+c---------------------------------
+ subroutine get_is_max(idum,ifar,n,imade,k)
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CSA'
+ include 'COMMON.BANK'
+ double precision emax
+
+ do i=imade+1,n
+ emax=-9.d190
+ do j=1,nbank
+ if(ibank(j).eq.k .and. bene(j).gt.emax) then
+ emax=bene(j)
+ is(i)=j
+ endif
+ enddo
+ call save_is(i)
+ enddo
+
+ return
+ end
projects / unres.git / blobdiff