C C $Id: resgrid.F,v 1.9 1998/07/16 16:40:38 jjv5 Exp arjan $ C C------------------------------------------------------------------------ subroutine resgrid(ipoint, natot, nbig, ibig, ierror) C C grid based subsetting C C written by Arjan van der Vaart, July-August 1997 C . January 1998 (pbc) C C C subsetting is residue-based for core and buffers. C implicit double precision (a-h,o-z) #include "divcon.dim" #include "divcon.h" C locals: logical nottaken dimension igsub1(maxatm), igsub2(maxatm), igsub3(maxatm), & igsub4(maxatm), nottaken(maxatm) CC-------------------------------------------------------------------CC C initialize ncell etc. This needs to be done everytime C this routine is called because it is dependent on the C maximum and minimum coordinates (can change during a NPT-MC-run) call setncell(ierror) if (ierror.ne.0) return C set up the grids: only residue-based is needed call mkrgrid C initialize some stuff before subsetting maxatsub = 0 imaxatsub = 0 n = 0 ngsub = 0 ngatot = 0 if (prtsub) write(iout,'(" SUBSYSTEM DEFINITIONS:")') C do the actual subsetting do 90 iz=0,ngstep(3) do 10 i=1,3 ibg(3,i) = iextl(3,i) + igstep(3)*iz ied(3,i) = ibg(3,i) + iextr(3,i) gridmin(3,i) = iz*gstep(3) + extl(3,i) gridmax(3,i) = gridmin(3,i) + extr(3,i) 10 continue do 80 iy=0,ngstep(2) do 20 i=1,3 ibg(2,i) = iextl(2,i) + igstep(2)*iy ied(2,i) = ibg(2,i) + iextr(2,i) gridmin(2,i) = iy*gstep(2) + extl(2,i) gridmax(2,i) = gridmin(2,i) + extr(2,i) 20 continue do 70 ix=0,ngstep(1) do 30 i=1,3 ibg(1,i) = iextl(1,i) + igstep(1)*ix ied(1,i) = ibg(1,i) + iextr(1,i) gridmin(1,i) = ix*gstep(1) + extl(1,i) gridmax(1,i) = gridmin(1,i) + extr(1,i) 30 continue k1 = 0 k2 = 0 k3 = 0 k4 = 0 C . core is residue-based, based on *real* distances call getrcorners(1,k1,igsub1,k2,igsub2) call getrribs(1,k1,igsub1,k2,igsub2) call getrplanes(1,k1,igsub1,k2,igsub2) call getrcore(k1,igsub1) if (combsub.and.(k1.ne.0)) then call cmbcore(k1,igsub1,k2,igsub2) endif if (k1.ne.0) then if (ldgbuff1) then C . buffers are residue-based if (pbc) then do i=1,natoms nottaken(i) = .true. enddo do i=1,k1 nottaken(igsub1(i)) = .false. enddo call pgetrcorners(2, k2, igsub2, k3, igsub3) call pgetrribs(2, k2, igsub2, k3, igsub3) call pgetrplanes(2, k2, igsub2, k3, igsub3) call pgetrinner(2, k2, igsub2) call busort(k2,igsub2,nottaken) else call getrcorners(2, k2, igsub2, k3, igsub3) call getrribs(2, k2, igsub2, k3, igsub3) call getrplanes(2, k2, igsub2, k3, igsub3) call getrinner(2, k2, igsub2) endif if (ldgbuff2) then if (pbc) then do i=1,k2 nottaken(igsub2(i)) = .false. enddo call pgetrcorners(3, k3, igsub3, k4, igsub4) call pgetrribs(3, k3, igsub3, k4, igsub4) call pgetrplanes(3, k3, igsub3, k4, igsub4) call pgetrinner(3, k3, igsub3) call busort(k3,igsub3,nottaken) else call getrcorners(3, k3, igsub3, k4, igsub4) call getrribs(3, k3, igsub3, k4, igsub4) call getrplanes(3, k3, igsub3, k4, igsub4) call getrinner(3, k3, igsub3) endif else k3 = 0 endif else k2 = 0 endif nsub = nsub + 1 ngsub = ngsub + 1 if (nsub .gt. maxsub) then write(iout,'(" MAXIMUM NUMBER OF ", & "SUBSYSTEMS (",I5,") REACHED." & /" CHANGE MAXSUB IN divcon.dim AND ", & "RECOMPILE.")') maxsub ierror = 1 return endif n = k1 + k2 + k3 natot = natot + n ngatot = ngatot + n if (natot .gt. mslist) then write(iout,'(/" ATOM LIST STORAGE LIMIT ", & "REACHED.",/" INCREASE MSLIST PARAMETER ", & "IN divcon.dim AND RECOMPILE")') ierror = 1 return endif if (n.gt.maxatsub) then maxatsub = n imaxatsub = nsub if (n.gt.nbig) then nbig = n ibig = nsub endif endif if (prtsub) call printsub(nsub, k1, igsub1, & k2, igsub2, k3, igsub3) iatom1(nsub) = ipoint do 40 iat=1,k1 iatoms(ipoint) = igsub1(iat) iabuff(ipoint) = 0 ipoint = ipoint + 1 40 continue do 50 iat=1,k2 iatoms(ipoint) = igsub2(iat) iabuff(ipoint) = 1 ipoint = ipoint + 1 50 continue do 60 iat=1,k3 iatoms(ipoint) = igsub3(iat) iabuff(ipoint) = 2 ipoint = ipoint + 1 60 continue C . sort atom list and carry along buffer status istart = iatom1(nsub) call bsort(n,iatoms(istart),iabuff(istart)) endif 70 continue 80 continue 90 continue iatom1(nsub+1) = ipoint natavg = ngatot/ngsub C-RDC write(iout,'(//" SUMMARY OF RESIDUE WISE, GRID BASED SUBSETTING:", C-RDC & /" ----------------------------------------------", C-RDC & //" EFFECTIVE BOXDIMENSIONS FOR SUBSETTING:", C-RDC & /" X = ",F10.5," (",F10.5," - ",F10.5,")", C-RDC & /" Y = ",F10.5," (",F10.5," - ",F10.5,")", C-RDC & /" Z = ",F10.5," (",F10.5," - ",F10.5,")")') C-RDC & gdxyz(1),gxyzmin(1),gxyzmax(1), C-RDC & gdxyz(2),gxyzmin(2),gxyzmax(2), C-RDC & gdxyz(3),gxyzmin(3),gxyzmax(3) C-RDC write(iout,'(/" NUMBER OF SUBSYSTEMS = ",I4, C-RDC & /" AVERAGE SUBSYSTEM SIZE = ",I4, C-RDC & /" MAXIMUM SUBSYSTEM SIZE = ",I4, C-RDC & /" LARGEST SUBSYSTEM = #",I4/)') C-RDC & ngsub,natavg,maxatsub,imaxatsub if (screen) then C-RDC write(iscr,'(/"resgrid subsetting: nsub = ",i4,", av. size = ", C-RDC & i4," max. size = ",i4)') ngsub,natavg,maxatsub endif C in case eigenvectors will be written if (prtvec) then ncores = 1 isubend(0) = 0 isubend(1) = nsub endif end