Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update QCG input file processing #243

Merged
merged 1 commit into from
Dec 7, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 1 addition & 5 deletions src/classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -349,7 +349,7 @@ module crest_data
character(len=512) :: homedir !> original directory from which calculation was started
character(len=512) :: scratchdir !> path to the scratch directory
character(len=:),allocatable :: cmd
character(len=:),allocatable :: inputcoords,inputcoords_solv,inputcoords_solu
character(len=:),allocatable :: inputcoords
character(len=:),allocatable :: wbofile
character(len=:),allocatable :: atlist
character(len=:),allocatable :: chargesfilename
Expand Down Expand Up @@ -385,8 +385,6 @@ module crest_data

!>--- property data objects
type(protobj) :: ptb
type(protobj) :: ptb_solvent
type(protobj) :: ptb_solute

!>--- saved constraints
type(constra) :: cts
Expand All @@ -408,8 +406,6 @@ module crest_data

!>--- reference structure data (the input structure)
type(refdata) :: ref
type(refdata) :: qcg_solvent
type(refdata) :: qcg_solute

!>--- QCG data
integer :: qcg_runtype = 0 !> Default is grow, 1= ensemble & opt, 2= e_solv, 3= g_solv
Expand Down
171 changes: 7 additions & 164 deletions src/confparse.f90
Original file line number Diff line number Diff line change
Expand Up @@ -491,6 +491,12 @@ subroutine parseflags(env,arg,nra)

case ('-solvtool','-qcg')
! env%mddumpxyz = 1000
! Set solute file if present
if(i == 2) env%solu_file = trim(arg(i-1))
! Set solvent file if prensent
! If it is another argument, it doesent matter as solvent file is checke in solvtool
if (nra >= i+1) env%solv_file = trim(arg(i+1))
! Set QCG defaults
env%preopt = .false.
env%crestver = crest_solv
env%QCG = .true.
Expand All @@ -501,9 +507,6 @@ subroutine parseflags(env,arg,nra)
env%ewin = 3.0d0
env%doOHflip = .false. !Switch off OH-flip
if (env%iterativeV2) env%iterativeV2 = .false.
if (nra >= i+1) then
env%solv_file = trim(arg(i+1))
end if
exit

case ('-compress')
Expand Down Expand Up @@ -2369,8 +2372,7 @@ subroutine inputcoords(env,arg)

!>--- Redirect for QCG input reading
if (env%QCG) then
arg2 = env%solv_file
call inputcoords_qcg(env,arg,arg2)
!Input coordinates are processed in solvtool.f90 file during solvtool subroutine
return
end if
!>---
Expand Down Expand Up @@ -2444,162 +2446,3 @@ subroutine inputcoords(env,arg)

return
end subroutine inputcoords

!========================================================================================!
!> Convert given QCG coordinate files into (TM format)
!> If the input is in SDF format, document the info to convert the
!> final ensemble back into this format
!========================================================================================!
subroutine inputcoords_qcg(env,arg1,arg2)
use iso_fortran_env,only:wp => real64
use crest_data
use strucrd
use zdata
use iomod
implicit none

type(systemdata) :: env
character(len=*) :: arg1,arg2

logical :: ex11,ex12,ex21,ex22,solu,solv
character(len=:),allocatable :: inputfile
type(coord) :: mol
type(zmolecule) :: zmol,zmol1
integer :: i

!--------------------Checking for input-------------!

inquire (file='solute',exist=solu)
inquire (file=arg1,exist=ex11)
inquire (file='coord',exist=ex12)

inquire (file='solvent',exist=solv)
inquire (file=arg2,exist=ex21)
if (len_trim(arg2) .eq. 0) then !Check if the second argument is just empty
ex21 = .false.
end if
inquire (file='coord',exist=ex22)

if (.not.ex11.and..not.ex12.and..not.solu) then
error stop 'No (valid) solute file! exit.'
else if (.not.ex21.and..not.ex22.and..not.solv) then
error stop 'No (valid) solvent file! exit.'
end if

!---------------Handling solute---------------------!

if (ex11.and.arg1(1:1) .ne. '-') then
call mol%open(arg1)
call mol%write('solute')
call mol%write('solute.xyz')
call mol%deallocate()
inputfile = trim(arg1)
write (*,*) 'Solute-file: ',arg1
env%solu_file = arg1
else if (solu.and..not.ex11) then
call copy('solute','solute.old')
inputfile = 'solute'
write (*,'(/,1x,a)') 'Solute-file: solute'
env%solu_file = 'solute'
else if (ex12.and..not.ex11) then !-- save coord as reference
call copy('coord','coord.old')
call copy('coord','solute')
write (*,'(/,1x,a)') 'Solute-file: coord'
inputfile = 'coord'
env%solu_file = 'coord'
else
write (*,*) 'An error occured processing the solute-input file'
end if

!--- if the input was a SDF file, special handling
env%sdfformat = .false.
call checkcoordtype(inputfile,i)
if (i == 31.or.i == 32) then
! call inpsdf(env,inputfile)
write (*,*) 'QCG currently does not support SDF-files.'
end if

!---------------Handling solvent---------------------!

if (ex21.and.arg2(1:1) .ne. '-'.and.arg2(1:1) .ne. ' ') then
call mol%open(arg2)
call mol%write('solvent')
call mol%write('solvent.xyz')
call mol%deallocate()
inputfile = trim(arg2)
write (*,*) 'Solvent-file: ',arg2
env%solv_file = arg2
else if (solv.and..not.ex21) then
call copy('solvent','solvent.old')
inputfile = 'solvent'
write (*,'(/,1x,a)') 'Solvent-file: solvent'
env%solv_file = 'solvent'
else if (ex22.and..not.ex21) then !-- save coord as reference
call copy('coord','coord.old')
call copy('coord','solvent')
inputfile = 'coord'
write (*,'(/,1x,a)') 'Solvent-file: coord'
env%solv_file = 'coord'
else
write (*,*) 'An error occured during solvent-file processing'
end if

!--- if the input was a SDF file, special handling
env%sdfformat = .false.
call checkcoordtype(inputfile,i)
if (i == 31.or.i == 32) then
!call inpsdf(env,inputfile)
write (*,*) 'QCG currently does not support SDF-files.'
end if

!-------------Checking both files and saving them---------------------------!

!--- after this point there should always be a solvent and solute file present
if (.not.allocated(env%inputcoords_solu)) env%inputcoords_solu = 'solute'
if (.not.allocated(env%inputcoords_solv)) env%inputcoords_solv = 'solvent'

call mol%open('solute')

env%nat = mol%nat
!--- solute geo
env%qcg_solute%nat = mol%nat
env%qcg_solute%at = mol%at
env%qcg_solute%xyz = mol%xyz
if (any((/crest_mfmdgc,crest_imtd,crest_imtd2/) == env%crestver)) then
if (.not.env%autozsort) then
env%qcg_solute%ntopo = mol%nat*(mol%nat+1)/2
allocate (env%qcg_solute%topo(env%qcg_solute%ntopo))
call quicktopo(mol%nat,mol%at,mol%xyz,env%qcg_solute%ntopo,env%qcg_solute%topo)
end if
end if
call mol%deallocate()

call mol%open('solvent')

env%nat = mol%nat
env%rednat = env%nat !get the number of atoms and the reduced number of atoms if some of them are excluded from the RMSD calc in V2
!--- solvent geo
env%qcg_solvent%nat = mol%nat
env%qcg_solvent%at = mol%at
env%qcg_solvent%xyz = mol%xyz
if (any((/crest_mfmdgc,crest_imtd,crest_imtd2/) == env%crestver)) then
if (.not.env%autozsort) then
env%qcg_solvent%ntopo = mol%nat*(mol%nat+1)/2
allocate (env%qcg_solvent%topo(env%qcg_solvent%ntopo))
call quicktopo(mol%nat,mol%at,mol%xyz,env%qcg_solvent%ntopo,env%qcg_solvent%topo)
end if
end if
call mol%deallocate()

!--- for protonation/deprotonation applications get ref. number of fragments
!--- also get some other structure based info
call simpletopo_file('solute',zmol,.false.,.false.,'')
env%ptb_solute%nfrag = zmol%nfrag
call zmol%deallocate()

call simpletopo_file('solvent',zmol1,.false.,.false.,'')
env%ptb_solvent%nfrag = zmol1%nfrag
call zmol1%deallocate()

return
end subroutine inputcoords_qcg
10 changes: 1 addition & 9 deletions src/crest_main.f90
Original file line number Diff line number Diff line change
Expand Up @@ -88,15 +88,7 @@ program CREST
!> SOME I/O STUFF
!=========================================================================================!
!>--- check for the coord file in the working directory
if (env%crestver .eq. crest_solv) then
inquire (file='solute',exist=ex1)
inquire (file='solvent',exist=ex2)
if (.not.ex1) then
error stop 'No solute file found. Exit.'
else if (.not.ex2) then
error stop 'No solvent file found. Exit.'
end if
else
if (env%crestver /= crest_solv) then
inquire (file='coord',exist=ex)
if (.not.ex) then
error stop 'No coord file found. Exit.'
Expand Down
Loading
Loading