Update GFN-FF submodule (crest-lab#266)

* Runtype for performing singlepoints for an entire ensemble
* Update gfnff git submodule

---------

Signed-off-by: Philipp Pracht <pp555@cam.ac.uk>

src/algos/ConfSolv.F90

@@ -99,8 +99,8 @@ subroutine cs_shutdown(io)
       write (stdout,'(/,a,i0)') 'Shutting down http://localhost/',cs_port
       call kill(cs_pid,9,io)
       deallocate (cs_pid)
+      call cs_shutdown2(io)
     end if
-    call cs_shutdown2(io)
 
   end subroutine cs_shutdown
 

src/algos/singlepoint.f90

@@ -230,3 +230,96 @@ subroutine crest_xtbsp(env,xtblevel,molin)
   call mol%deallocate()
 end subroutine crest_xtbsp
 
+!========================================================================================!
+!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<!
+!========================================================================================!
+subroutine crest_ensemble_singlepoints(env,tim)
+!***********************************************
+!* subroutine crest_ensemble_singlepoints
+!* This routine implements a standalone runtype
+!* to perform singlepoint evaluations along an
+!* ensemble or trajectory file.
+!***********************************************
+  use crest_parameters,only:wp,stdout,bohr,angstrom
+  use crest_data
+  use crest_calculator
+  use strucrd
+  use optimize_module
+  use utilities, only: dumpenergies
+  implicit none
+  type(systemdata),intent(inout) :: env
+  type(timer),intent(inout)      :: tim
+  type(coord) :: mol,molnew
+  integer :: i,j,k,l,io,ich,c
+  logical :: pr,wr,ex
+!========================================================================================!
+  type(calcdata) :: calc
+
+  real(wp) :: energy,gnorm
+  real(wp),allocatable :: grad(:,:)
+
+  character(len=:),allocatable :: ensnam
+  integer :: nat,nall
+  real(wp),allocatable :: eread(:)
+  real(wp),allocatable :: xyz(:,:,:)
+  integer,allocatable  :: at(:)
+  character(len=80) :: atmp
+  real(wp) :: percent
+  character(len=52) :: bar
+!========================================================================================!
+ write (*,*)
+!>--- check for the ensemble file
+  inquire (file=env%ensemblename,exist=ex)
+  if (ex) then
+    ensnam = env%ensemblename
+  else
+    write (stdout,*) 'no ensemble file provided.'
+    return
+  end if
+
+!>--- start the timer
+  call tim%start(14,'Ensemble singlepoints')
+
+!>---- read the input ensemble
+  call rdensembleparam(ensnam,nat,nall)
+  if (nall .lt. 1) return
+  allocate (xyz(3,nat,nall),at(nat),eread(nall))
+  call rdensemble(ensnam,nat,nall,at,xyz,eread)
+!>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<!
+!>--- Important: crest_oloop requires coordinates in Bohrs
+  xyz = xyz / bohr
+!>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<!
+
+!>--- set OMP parallelization
+  if (env%autothreads) then
+    !>--- usually, one thread per xtb job
+    call ompautoset(env%threads,7,env%omp,env%MAXRUN,nall)
+  end if
+!========================================================================================!
+  !>--- printout header
+  write (stdout,*)
+  write (stdout,'(10x,"┍",49("━"),"┑")')
+  write (stdout,'(10x,"│",14x,a,14x,"│")') "ENSEMBLE SINGLEPOINTS"
+  write (stdout,'(10x,"┕",49("━"),"┙")')
+  write (stdout,*)
+  write (stdout,'(1x,a,i0,a,1x,a)') 'Evaluationg all ',nall,' structures of file',trim(ensnam)
+  !>--- call the loop
+  call crest_sploop(env,nat,nall,at,xyz,eread,.true.)
+
+!>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<!
+!>--- Important: ensemble file must be written in AA
+  xyz = xyz / angstrom
+!>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<!
+!>--- write output ensemble
+  call wrensemble(ensemblefile,nat,nall,at,xyz,eread)
+  write(stdout,'(/,a,a,a)') 'Ensemble with updated energies written to <',ensemblefile,'>'
+
+  call dumpenergies('crest.energies',eread)
+  write(stdout,'(/,a,a,a)') 'List of energies written to <','crest.energies','>' 
+
+  deallocate (eread,at,xyz)
+!========================================================================================!
+  call tim%stop(14)
+  return
+end subroutine crest_ensemble_singlepoints
+

src/classes.f90

@@ -81,6 +81,7 @@ module crest_data
   integer,parameter,public :: crest_scanning   = 270
   integer,parameter,public :: crest_rigcon     = 271
   integer,parameter,public :: crest_trialopt   = 272
+  integer,parameter,public :: crest_ensemblesp = 273
 !>> <<!
   integer,parameter,public :: crest_test       = 456
 

src/confparse.f90

@@ -398,6 +398,20 @@ subroutine parseflags(env,arg,nra)
         env%inputcoords = env%ensemblename !just for a printout
         exit
 
+      case ('-mdsp','-ensemblesp') !> Singlepoints along ensemble
+        env%crestver = crest_ensemblesp
+        atmp = ''
+        env%preopt = .false.
+        env%ensemblename = 'none selected'
+        if (nra .ge. (i+1)) atmp = adjustl(arg(i+1))
+        if ((atmp(1:1) /= '-').and.(len_trim(atmp) .ge. 1)) then
+          env%ensemblename = trim(atmp)
+        end if
+        call xyz2coord(env%ensemblename,'coord') !> write coord from lowest structure
+        env%inputcoords = env%ensemblename !> just for a printout
+        exit
+
+
       case ('-pka','-pKa')  !> pKa calculation script
         env%crestver = crest_pka
         env%runver = 33
@@ -1037,7 +1051,9 @@ subroutine parseflags(env,arg,nra)
           env%mdstep = 1.5d0
           env%hmass = 5.0d0
           ctype = 5 !> bond constraint activated
-          bondconst = .true.
+          if (any((/crest_imtd,crest_imtd2/) == env%crestver)) then
+            bondconst = .true.
+          endif
           env%cts%cbonds_md = .true.
           env%checkiso = .true.
         case ('stereoisomers')

src/crest_main.f90

@@ -286,7 +286,10 @@ program CREST
     call crest_rigidconf(env,tim)
 
   case (crest_trialopt) !> test optimization standalone
-    call trialOPT(env)    
+    call trialOPT(env)
+
+  case (crest_ensemblesp) !> singlepoints along ensemble
+    call crest_ensemble_singlepoints(env,tim)    
 
   case (crest_test)
     call crest_playground(env,tim)

src/parsing/parse_maindata.f90

@@ -129,6 +129,10 @@ subroutine parse_main_c(env,key,val)
       case ('screen_ensemble','screen')
         env%preopt = .false.
         env%crestver = crest_screen
+      case ('ensemble_singlepoints','ensemblesp','mdsp')
+        env%preopt = .false.
+        env%crestver = crest_ensemblesp
+
       case ('md','mtd','metadynamics','dynamics')
         env%preopt = .false.
         env%crestver = crest_moldyn

src/utilmod.f90

@@ -37,6 +37,7 @@ module utilities
   public :: revlin
   public :: TRJappendto_skipfirst
   public :: XYZappendto
+  public :: dumpenergies
 
   !> functions
   public :: lin
@@ -519,6 +520,25 @@ subroutine XYZappendto(from,to)
     close (tunit)
   end subroutine XYZappendto
 
+!============================================================!
+
+  subroutine dumpenergies(filename,eread)
+!****************************
+!* write energies to a file
+!****************************
+    implicit none
+    character(len=*),intent(in) :: filename
+    real(wp),intent(in) :: eread(:)
+    integer :: ich,io,l,i
+
+    open (newunit=ich,file=filename)
+    l = size(eread,1)
+    do i = 1,l
+      write (ich,'(f25.15)') eread(i)
+    end do
+    close (ich)
+  end subroutine dumpenergies
+
 !========================================================================================!
 !========================================================================================!
 end module utilities