From a335e5ed2b6e77aef321d20ef05e137f8c1e2c5c Mon Sep 17 00:00:00 2001 From: Philipp Pracht <42216327+pprcht@users.noreply.github.com> Date: Fri, 16 Feb 2024 14:27:26 +0000 Subject: [PATCH] Update GFN-FF submodule (#266) * Runtype for performing singlepoints for an entire ensemble * Update gfnff git submodule --------- Signed-off-by: Philipp Pracht --- src/algos/ConfSolv.F90 | 2 +- src/algos/singlepoint.f90 | 93 ++++++++++++++++++++++++++++++++++ src/classes.f90 | 1 + src/confparse.f90 | 18 ++++++- src/crest_main.f90 | 5 +- src/parsing/parse_maindata.f90 | 4 ++ src/utilmod.f90 | 20 ++++++++ subprojects/gfnff | 2 +- 8 files changed, 141 insertions(+), 4 deletions(-) diff --git a/src/algos/ConfSolv.F90 b/src/algos/ConfSolv.F90 index 29ae231e..91a8dd66 100644 --- a/src/algos/ConfSolv.F90 +++ b/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 diff --git a/src/algos/singlepoint.f90 b/src/algos/singlepoint.f90 index 94375e71..bba47042 100644 --- a/src/algos/singlepoint.f90 +++ b/src/algos/singlepoint.f90 @@ -230,3 +230,96 @@ subroutine crest_xtbsp(env,xtblevel,molin) call mol%deallocate() end subroutine crest_xtbsp +!========================================================================================! +!>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<--- 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 + diff --git a/src/classes.f90 b/src/classes.f90 index 254c4182..800bc966 100644 --- a/src/classes.f90 +++ b/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 !>> < 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') diff --git a/src/crest_main.f90 b/src/crest_main.f90 index 5ce51a6e..108404a4 100644 --- a/src/crest_main.f90 +++ b/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) diff --git a/src/parsing/parse_maindata.f90 b/src/parsing/parse_maindata.f90 index 219d72d0..bad81ed5 100644 --- a/src/parsing/parse_maindata.f90 +++ b/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 diff --git a/src/utilmod.f90 b/src/utilmod.f90 index ac8402a7..73836b0b 100644 --- a/src/utilmod.f90 +++ b/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 diff --git a/subprojects/gfnff b/subprojects/gfnff index 9ea8abdc..4d68ac1a 160000 --- a/subprojects/gfnff +++ b/subprojects/gfnff @@ -1 +1 @@ -Subproject commit 9ea8abdcad146830bbff967f069532b15dfe7af4 +Subproject commit 4d68ac1ab8df5999d3493715a86c13786bac6bfb