--- /dev/null
+ module map_
+!-----------------------------------------------------------------------------
+ use io_units
+ use names
+ use geometry_data
+ use control_data
+ use energy_data
+ use map_data
+ implicit none
+!-----------------------------------------------------------------------------
+!
+!
+!-----------------------------------------------------------------------------
+ contains
+!-----------------------------------------------------------------------------
+! check_sc_map.f
+!-----------------------------------------------------------------------------
+ subroutine check_sc_map
+! Subroutine is checking if the fitted function which describs sc_rot_pot
+! is correct, printing, alpha,beta, energy, data - for some known theta.
+! theta angle is read from the input file. Sc_rot_pot are printed
+! for the second residue in sequance.
+ use geometry, only:chainbuild
+ use energy, only:vec_and_deriv,esc
+! include 'DIMENSIONS'
+! include 'COMMON.VAR'
+! include 'COMMON.GEO'
+! include 'COMMON.INTERACT'
+ real(kind=8) :: xx,yy,zz,al,om
+ real(kind=8) :: escloc, escloc_min
+ real(kind=8),dimension(50000) :: escloc_ene,alph_plot,beta_plot
+ integer,dimension(5000) :: al_plot,be_plot
+ integer :: iialph, iibet,it
+
+!el local variables
+ integer :: i,j
+
+ write (2,*) "Side-chain-rotamer potential energy map!!!!"
+ escloc_min = 1000000.00
+! it=itype(2)
+ i = 0
+ do iialph=0,18
+ do iibet=-18,18
+ i = i + 1
+ al = iialph*10.0d0*deg2rad
+ om = iibet*10.0d0*deg2rad
+ zz = dcos(al)
+ xx = -dsin(al)*dcos(om)
+ yy = -dsin(al)*dsin(om)
+ alph(2)=dacos(xx)
+ omeg(2)=-datan2(zz,yy)
+ al_plot(i)=alph(2)*rad2deg
+ be_plot(i)=omeg(2)*rad2deg
+! write(2,*) alph(2)*rad2deg, omeg(2)*rad2deg
+ alph_plot(i) = al*rad2deg
+ beta_plot(i) = om*rad2deg
+ call chainbuild
+ call vec_and_deriv
+ call esc(escloc)
+ escloc_ene(i) = escloc
+ if (escloc_min.gt.escloc_ene(i)) escloc_min=escloc_ene(i)
+ enddo
+ enddo
+! write (2,*) "escloc_min = ", escloc_min
+ print *,"i",i
+ do j = 1,i
+ write (2,'(3f10.3,2i9,f12.5)') alph_plot(j), &
+ beta_plot(j),theta(3)*rad2deg, al_plot(j),be_plot(j),&
+ escloc_ene(j) !- escloc_min
+ enddo
+ return
+ end subroutine check_sc_map
+!-----------------------------------------------------------------------------
+! map.f
+!-----------------------------------------------------------------------------
+ subroutine map
+
+ use geometry, only:chainbuild,geom_to_var
+ use energy
+ use minimm, only:minimize
+! implicit real*8 (a-h,o-z)
+! include 'DIMENSIONS'
+! include 'COMMON.MAP'
+! include 'COMMON.VAR'
+! include 'COMMON.GEO'
+! include 'COMMON.DERIV'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.NAMES'
+! include 'COMMON.CONTROL'
+! include 'COMMON.TORCNSTR'
+ real(kind=8) :: energia(0:n_ene)
+ character(len=5) :: angid(4)=reshape((/'PHI ','THETA','ALPHA',&
+ 'OMEGA'/),shape(angid))
+ real(kind=8) :: ang_list(10)
+!el real(kind=8),dimension(6*nres) :: g,x !(maxvar) (maxvar=6*maxres)
+ real(kind=8),dimension(:),allocatable :: g,x !(maxvar) (maxvar=6*maxres)
+ integer :: nn(10)
+!el local variables
+ integer :: i,iii,ii,j,k,nmax,ntot,nf,iretcode,nfun
+ real(kind=8) :: etot,gnorm!,fdum
+ integer,dimension(1) :: uiparm
+ real(kind=8),dimension(1) :: urparm
+ allocate(x(6*nres),g(6*nres))
+
+ write (iout,'(a,i3,a)')'Energy map constructed in the following ',&
+ nmap,' groups of variables:'
+ do i=1,nmap
+ write (iout,'(2a,i3,a,i3)') angid(kang(i)),' of residues ',&
+ res1(i),' to ',res2(i)
+ enddo
+ nmax=nstep(1)
+ do i=2,nmap
+ if (nmax.lt.nstep(i)) nmax=nstep(i)
+ enddo
+ ntot=nmax**nmap
+ iii=0
+ write (istat,'(1h#,a14,29a15)') (" ",k=1,nmap),&
+ (ename(print_order(k)),k=1,nprint_ene),"ETOT","GNORM"
+ do i=0,ntot-1
+ ii=i
+ do j=1,nmap
+ nn(j)=mod(ii,nmax)+1
+ ii=ii/nmax
+ enddo
+ do j=1,nmap
+ if (nn(j).gt.nstep(j)) goto 10
+ enddo
+ iii=iii+1
+!d write (iout,*) i,iii,(nn(j),j=1,nmap)
+ do j=1,nmap
+ ang_list(j)=ang_from(j) &
+ +(nn(j)-1)*(ang_to(j)-ang_from(j))/nstep(j)
+ do k=res1(j),res2(j)
+ goto (1,2,3,4), kang(j)
+ 1 phi(k)=deg2rad*ang_list(j)
+ if (minim) phi0(k-res1(j)+1)=deg2rad*ang_list(j)
+ goto 5
+ 2 theta(k)=deg2rad*ang_list(j)
+ goto 5
+ 3 alph(k)=deg2rad*ang_list(j)
+ goto 5
+ 4 omeg(k)=deg2rad*ang_list(j)
+ 5 continue
+ enddo ! k
+ enddo ! j
+ call chainbuild
+ if (minim) then
+ call geom_to_var(nvar,x)
+ call minimize(etot,x,iretcode,nfun)
+ print *,'SUMSL return code is',iretcode,' eval ',nfun
+! call intout
+ else
+ call zerograd
+ call geom_to_var(nvar,x)
+ endif
+ call etotal(energia)
+ etot = energia(0)
+ nf=1
+ nfl=3
+ call gradient(nvar,x,nf,g,uiparm,urparm,fdum)
+ gnorm=0.0d0
+ do k=1,nvar
+ gnorm=gnorm+g(k)**2
+ enddo
+ etot=energia(0)
+
+ gnorm=dsqrt(gnorm)
+! write (iout,'(6(1pe15.5))') (ang_list(k),k=1,nmap),etot,gnorm
+ write (istat,'(30e15.5)') (ang_list(k),k=1,nmap),&
+ (energia(print_order(ii)),ii=1,nprint_ene),etot,gnorm
+! write (iout,*) 'POINT',I,' ANGLES:',(ang_list(k),k=1,nmap)
+! call intout
+! call enerprint(energia)
+ 10 continue
+ enddo ! i
+! deallocate(x,g)
+ return
+ end subroutine map
+!-----------------------------------------------------------------------------
+ subroutine alloc_map_arrays
+
+! commom.map
+! common /mapp/
+ allocate(kang(nmap),res1(nmap),res2(nmap),nstep(nmap)) !(maxvar)
+ allocate(ang_from(nmap),ang_to(nmap)) !(maxvar)
+
+ return
+ end subroutine alloc_map_arrays
+!-----------------------------------------------------------------------------
+!-----------------------------------------------------------------------------
+ end module map_
projects / unres4.git / blobdiff