Add "--tmcosmo" mode for writing .cosmo files with TM convention.

src/prog/main.F90

@@ -1595,6 +1595,25 @@ subroutine parseArguments(env, args, inputFile, paramFile, accuracy, lgrad, &
             call env%error("No solvent name provided for COSMO", source)
          end if
 
+      case('--tmcosmo')
+         call args%nextArg(sec)
+         if (allocated(sec)) then
+            call set_gbsa(env, 'solvent', sec)
+            call set_gbsa(env, 'tmcosmo', 'true')
+            call args%nextArg(sec)
+            if (allocated(sec)) then
+               if (sec == 'reference') then
+                  gsolvstate = 1
+               else if (sec == 'bar1M') then
+                  gsolvstate = 2
+               else
+                  call env%warning("Unknown reference state '"//sec//"'", source)
+               end if
+            end if
+         else
+            call env%error("No solvent name provided for COSMO", source)
+         end if
+
       case('--cpcmx')
          if (get_xtb_feature('cpcmx')) then
             call args%nextArg(sec)

src/scf_module.F90

@@ -854,7 +854,7 @@ subroutine scf(env, mol, wfn, basis, pcem, xtbData, solvation, &
          type is (TCosmo)
             call open_file(ich, "xtb.cosmo", 'w')
             call solvation%writeCosmoFile(ich, mol%at, mol%sym, mol%xyz, &
-               & wfn%q, eel + ep + exb + merge(embd, ed + embd, allocated(scD4)))
+               & wfn%q, eel + ep + exb + merge(embd, ed + embd, allocated(scD4)),solvation%tmcosmo)
             call close_file(ich)
          end select
       end if

src/set_module.f90

@@ -2020,6 +2020,7 @@ subroutine set_gbsa(env,key,val)
    logical,save :: set6 = .true.
    logical,save :: set7 = .true.
    logical,save :: set8 = .true.
+   logical,save :: set9 = .true.
    select case(key)
    case default ! do nothing
       call env%warning("the key '"//key//"' is not recognized by gbsa",source)
@@ -2073,8 +2074,15 @@ subroutine set_gbsa(env,key,val)
    case('cosmo')
       if (getValue(env,val,ldum).and.set7) set%solvInput%cosmo = ldum
       set7 = .false.
+   case('tmcosmo')
+      if (getValue(env,val,ldum).and.set8) then
+         set%solvInput%cosmo = ldum
+         set%solvInput%tmcosmo = .true.
+      end if
+      set8 = .false.
    case('cpcmx')
-      if (set8) set%solvInput%cpxsolvent = val
+      if (set9) set%solvInput%cpxsolvent = val
+      set9 = .false.
    end select
 end subroutine set_gbsa
 

src/solv/cosmo.f90

@@ -95,6 +95,9 @@ module xtb_solv_cosmo
       real(wp), allocatable :: dsdr(:, :)
       real(wp), allocatable :: dsdrt(:, :, :)
 
+      !> TM convention?
+      logical :: tmcosmo = .false.
+
    contains
 
       !> Update coordinates and internal state
@@ -826,7 +829,7 @@ subroutine update_nnlist_sasa(nat, xyz, srcut, nnsas, nnlists)
 end subroutine update_nnlist_sasa
 
 !> Write a COSMO file output
-subroutine writeCosmoFile(self, unit, num, sym, xyz, qat, energy)
+subroutine writeCosmoFile(self, unit, num, sym, xyz, qat, energy, tmcosmo)
 
    !> COSMO container
    class(TCosmo), intent(in) :: self
@@ -849,11 +852,20 @@ subroutine writeCosmoFile(self, unit, num, sym, xyz, qat, energy)
    !> Total energy
    real(wp), intent(in) :: energy
 
+   !> Switch to TM convention for cosmo file output
+   logical, intent(in), optional :: tmcosmo
+
    integer :: ii, ig, iat
    real(wp) :: dielEnergy, keps
    real(wp), allocatable :: phi(:), zeta(:), area(:)
+   logical :: tm
 
    allocate(phi(self%ddCosmo%ncav), zeta(self%ddCosmo%ncav), area(self%ddCosmo%ncav))
+   if (present(tmcosmo)) then
+      tm = tmcosmo
+   else
+      tm = .false.
+   end if
    ! Reset potential on the cavity, note that the potential is expected in e/Å
    call getPhi(qat, self%jmat, phi)
    ii = 0
@@ -870,6 +882,9 @@ subroutine writeCosmoFile(self, unit, num, sym, xyz, qat, energy)
       end do
    end do
 
+   !! Switch convention for TM mode
+   if (tm) zeta=-zeta
+
 
     ! Dielectric energy is the energy on the dielectric continuum
    keps = 0.5_wp * (1.0_wp - 1.0_wp/self%dielectricConst)
@@ -940,6 +955,7 @@ subroutine writeCosmoFile(self, unit, num, sym, xyz, qat, energy)
       & "#", &
       & "#"
 
+
    ii = 0
    do iat = 1, self%ddCosmo%nat
       do ig = 1, self%ddCosmo%ngrid

src/solv/input.F90

@@ -42,6 +42,9 @@ module xtb_solv_input
       !> Use COSMO solvation model instead of Born based one
       logical :: cosmo = .false.
 
+      !> Use TM convention for COSMO solvation model
+      logical :: tmcosmo = .false.
+
       !> Additional CPCM-X mode
       character(len=:),allocatable :: cpxsolvent
 

src/solv/model.f90

@@ -693,6 +693,7 @@ subroutine newSolvationModel(self, env, model, num)
       allocate(cosmo)
       call init_(cosmo, env, num, self%dielectricConst, self%nAng, self%bornScale, &
          & self%vdwRad, self%surfaceTension, self%probeRad, srcut)
+      cosmo%tmcosmo=self%tmcosmo
       call move_alloc(cosmo, model)
    else
       allocate(born)

src/xhelp.f90

@@ -165,6 +165,9 @@ subroutine help(iunit)
    "    Additionally, the dielectric constant can be set manually or an ideal conductor", &
    "    can be chosen by setting epsilon to infinity.",&
    "",&
+   "--tmcosmo SOLVENT/EPSILON",&
+   "    same as --cosmo, but uses TM convention for writing the .cosmo files.",&
+   "",&
    "--cpcmx SOLVENT",&
    "    extended conduction-like polarizable continuum solvation model (CPCM-X),",&
    "    available solvents are all solvents included in the Minnesota Solvation Database.",&