Conditional cosmo file generation.

This commit is contained in:
rpseng 2025-03-13 18:03:18 -03:00
parent 2c7770f24e
commit d61602b6cb
2 changed files with 48 additions and 34 deletions

View file

@ -79,6 +79,7 @@ c
double precision dielec,dielecinf,rsolv,rad,zeta,gammas,swtol
character*255 field
logical do_cosmo
logical do_cosmo_file
logical do_gasphase
logical status
character*30 tag
@ -326,6 +327,15 @@ c
$ do_cosmo_model))
$ call errquit('cosmo_input: rtdb put failed',911, RTDB_ERR)
c
c --- do_cosmo_file flag
c
else if(inp_compare(.false.,'do_cosmo_file',field)) then
do_cosmo_file=.true.
if (.not.
$ rtdb_put(rtdb,'cosmo:do_cosmo_file',mt_log,1,do_cosmo_file))
$ call errquit('cosmo_input: cannot put do_cosmo_file in rtdb',
$ 0,rtdb_err)
c
c --- rsolv: Used in Klamt-Schuurmann, but not in York-Karplus
c --- approach
c

View file

@ -479,7 +479,7 @@ c
double precision efcs(nefc) !< [Input] the COSMO charge
!< surface area for KS
double precision msrf(nefc) !< [Input] the COSMO charge
!< surface area for YS
!< surface area for YK
c
double precision efcz(nefc) !< [Output] the COSMO charges
double precision efczz(nefc) !< [Input[ the COSMO zeta value for
@ -490,6 +490,7 @@ c
double precision ecos !< [Output] the energy contribution due to
!< the COSMO charges
logical status
logical do_cosmo_file
c
logical dbug,more,out,direct,noall,all,elec,nucl,iefc_done
character*16 at_tag
@ -1219,52 +1220,55 @@ c
c
c ----- printing cosmo file for COSMO-RS/SAC codes -----
c
cosmo_file_in = "cosmo"
call util_file_name(cosmo_file_in,.false.,.false.,cosmo_file)
call util_file_name_resolve(cosmo_file,.false.)
if(ga_nodeid().eq.0) then
if(.not.util_io_unit(80,90,fn))
+ call errquit('cannot get free unit', 0,
+ 0)
if (.not. rtdb_get(rtdb,'cosmo:do_cosmo_file',mt_log,1,
$ do_cosmo_file)) do_cosmo_file=.false.
if (do_cosmo_file) then
cosmo_file_in = "cosmo"
call util_file_name(cosmo_file_in,.false.,.false.,cosmo_file)
call util_file_name_resolve(cosmo_file,.false.)
if(ga_nodeid().eq.0) then
if(.not.util_io_unit(80,90,fn))
+ call errquit('cannot get free unit', 0, 0)
c
open(unit=fn,form="formatted",file=cosmo_file)
if (dbug) then
write(*,*) "printing cosmo file for COSMO-RS/SAC",
open(unit=fn,form="formatted",file=cosmo_file)
if (dbug) then
write(*,*) "printing cosmo file for COSMO-RS/SAC",
+ cosmo_file
end if
write(fn,"(a)") "$cosmo"
write(fn,*) " epsilon= ", dielec
write(fn,*) " rsolv= ", rsolv
end if
write(fn,"(a)") "$cosmo"
write(fn,*) " epsilon= ", dielec
write(fn,*) " rsolv= ", rsolv
c
write(fn,"(a)") "$cosmo_data"
write(fn,*) " nps= ", nefc
write(fn,*) " area= ", srfmol/(cau2ang**2)
write(fn,*) " volume= ", volmol/(cau2ang**3)
write(fn,"(a)") "$cosmo_data"
write(fn,*) " nps= ", nefc
write(fn,*) " area= ", srfmol/(cau2ang**2)
write(fn,*) " volume= ", volmol/(cau2ang**3)
c
write(fn,"(a)") "$coord_rad"
write(fn,"(a)") "#atom x y z [a.u.] element radius [A]"
do i = 1,nat
write(fn,*) i,
write(fn,"(a)") "$coord_rad"
write(fn,"(a)") "#atom x y z [a.u.] element radius [A]"
do i = 1,nat
write(fn,*) i,
+ dbl_mb(k_xyzpt +3*(i-1)),
+ dbl_mb(k_xyzpt+1+3*(i-1)),
+ dbl_mb(k_xyzpt+2+3*(i-1)),
+ int(dbl_mb(k_zanpt+i-1)),
+ ratm(i)*cau2ang
enddo ! i
enddo ! i
c
write(fn,"(a)") "$segment_information"
write(fn,"(a)") "# n atom x y z charge area charge/area"
write(fn,"(a)") "# coordinates [a.u.], area [A2]"
do ief=1,nefc
efcsi = efcs(ief)*(cau2ang**2)
if (do_cosmo_model.eq.DO_COSMO_YK)
+ efcsi = msrf(ief)*(cau2ang**2)
write(fn,"(a)") "$segment_information"
write(fn,"(a)") "# n atom x y z charge area charge/area"
write(fn,"(a)") "# coordinates [a.u.], area [A2]"
do ief=1,nefc
efcsi = efcs(ief)*(cau2ang**2)
if (do_cosmo_model.eq.DO_COSMO_YK)
+ efcsi = msrf(ief)*(cau2ang**2)
write(fn,*) ief, efciat(ief),
write(fn,*) ief, efciat(ief),
+ efcc(1,ief), efcc(2,ief), efcc(3,ief),
+ efcz(ief), efcsi, efcz(ief)/efcsi
end do
close(fn)
end do
close(fn)
end if
end if
c
c MN solvation models -->