compute_thermochemistry Subroutine

public subroutine compute_thermochemistry(coords, atomic_numbers, frequencies, n_atoms, n_freqs, result, temperature, pressure, symmetry_number, spin_multiplicity)

Main driver for thermochemistry calculations.

Computes all thermodynamic quantities from molecular geometry and vibrational frequencies.

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in) :: coords(:,:)

Coordinates (3, n_atoms) in Bohr
integer, intent(in) :: atomic_numbers(:)

Atomic numbers
real(kind=dp), intent(in) :: frequencies(:)

Frequencies in cm^-1
integer, intent(in) :: n_atoms

Number of atoms
integer, intent(in) :: n_freqs

Number of frequencies
type(thermochemistry_result_t), intent(out) :: result

Output results
real(kind=dp), intent(in), optional :: temperature

Temperature in K (default 298.15)
real(kind=dp), intent(in), optional :: pressure

Pressure in atm (default 1.0)
integer, intent(in), optional :: symmetry_number

Symmetry number (default 1)
integer, intent(in), optional :: spin_multiplicity

Spin multiplicity (default 1)

Calls

proc~~compute_thermochemistry~~CallsGraph proc~compute_thermochemistry compute_thermochemistry proc~compute_electronic_entropy compute_electronic_entropy proc~compute_thermochemistry->proc~compute_electronic_entropy proc~compute_moments_of_inertia compute_moments_of_inertia proc~compute_thermochemistry->proc~compute_moments_of_inertia proc~compute_partition_functions compute_partition_functions proc~compute_thermochemistry->proc~compute_partition_functions proc~compute_rotational_constants compute_rotational_constants proc~compute_thermochemistry->proc~compute_rotational_constants proc~compute_rotational_thermo compute_rotational_thermo proc~compute_thermochemistry->proc~compute_rotational_thermo proc~compute_translational_thermo compute_translational_thermo proc~compute_thermochemistry->proc~compute_translational_thermo proc~compute_vibrational_thermo compute_vibrational_thermo proc~compute_thermochemistry->proc~compute_vibrational_thermo proc~compute_zpe compute_zpe proc~compute_thermochemistry->proc~compute_zpe pic_syev pic_syev proc~compute_moments_of_inertia->pic_syev proc~element_mass element_mass proc~compute_moments_of_inertia->proc~element_mass to_char to_char proc~compute_moments_of_inertia->to_char warning warning proc~compute_moments_of_inertia->warning proc~compute_zpe->to_char proc~compute_zpe->warning

Called by

proc~~compute_thermochemistry~~CalledByGraph proc~compute_thermochemistry compute_thermochemistry proc~compute_mbe compute_mbe proc~compute_mbe->proc~compute_thermochemistry proc~print_vibrational_analysis print_vibrational_analysis proc~compute_mbe->proc~print_vibrational_analysis proc~gmbe_pie_coordinator gmbe_pie_coordinator proc~gmbe_pie_coordinator->proc~compute_thermochemistry proc~gmbe_pie_coordinator->proc~print_vibrational_analysis proc~hessian_coordinator hessian_coordinator proc~hessian_coordinator->proc~compute_thermochemistry proc~hessian_coordinator->proc~print_vibrational_analysis proc~print_vibrational_analysis->proc~compute_thermochemistry proc~serial_gmbe_pie_processor serial_gmbe_pie_processor proc~serial_gmbe_pie_processor->proc~compute_thermochemistry proc~serial_gmbe_pie_processor->proc~print_vibrational_analysis proc~unfragmented_calculation unfragmented_calculation proc~unfragmented_calculation->proc~compute_thermochemistry proc~unfragmented_calculation->proc~print_vibrational_analysis interface~hessian_coordinator hessian_coordinator interface~hessian_coordinator->proc~hessian_coordinator interface~unfragmented_calculation unfragmented_calculation interface~unfragmented_calculation->proc~unfragmented_calculation proc~compute_gmbe compute_gmbe proc~compute_gmbe->proc~print_vibrational_analysis proc~global_coordinator global_coordinator proc~global_coordinator->proc~compute_mbe proc~gmbe_run_distributed gmbe_context_t%gmbe_run_distributed proc~gmbe_run_distributed->proc~gmbe_pie_coordinator proc~gmbe_run_serial gmbe_context_t%gmbe_run_serial proc~gmbe_run_serial->proc~serial_gmbe_pie_processor proc~serial_fragment_processor serial_fragment_processor proc~serial_fragment_processor->proc~compute_mbe interface~global_coordinator global_coordinator interface~global_coordinator->proc~global_coordinator interface~serial_fragment_processor serial_fragment_processor interface~serial_fragment_processor->proc~serial_fragment_processor proc~distributed_unfragmented_hessian distributed_unfragmented_hessian proc~distributed_unfragmented_hessian->interface~hessian_coordinator proc~run_unfragmented_calculation run_unfragmented_calculation proc~run_unfragmented_calculation->interface~unfragmented_calculation interface~distributed_unfragmented_hessian distributed_unfragmented_hessian proc~run_unfragmented_calculation->interface~distributed_unfragmented_hessian interface~distributed_unfragmented_hessian->proc~distributed_unfragmented_hessian proc~mbe_run_distributed mbe_context_t%mbe_run_distributed proc~mbe_run_distributed->interface~global_coordinator proc~mbe_run_serial mbe_context_t%mbe_run_serial proc~mbe_run_serial->interface~serial_fragment_processor proc~run_calculation run_calculation proc~run_calculation->proc~run_unfragmented_calculation proc~compute_energy_and_forces compute_energy_and_forces proc~compute_energy_and_forces->proc~run_calculation proc~run_multi_molecule_calculations run_multi_molecule_calculations proc~run_multi_molecule_calculations->proc~run_calculation program~main main program~main->proc~run_calculation

Variables

Type Visibility Attributes Name Initial
real(kind=dp), private :: P
real(kind=dp), private :: S_total
real(kind=dp), private :: T
real(kind=dp), private :: center_of_mass(3)
integer, private :: mult
real(kind=dp), private :: principal_axes(3,3)
integer, private :: sigma

Source Code

   subroutine compute_thermochemistry(coords, atomic_numbers, frequencies, n_atoms, n_freqs, &
                                      result, temperature, pressure, symmetry_number, spin_multiplicity)
      !! Main driver for thermochemistry calculations.
      !!
      !! Computes all thermodynamic quantities from molecular geometry and vibrational frequencies.
      real(dp), intent(in) :: coords(:, :)           !! Coordinates (3, n_atoms) in Bohr
      integer, intent(in) :: atomic_numbers(:)       !! Atomic numbers
      real(dp), intent(in) :: frequencies(:)         !! Frequencies in cm^-1
      integer, intent(in) :: n_atoms                 !! Number of atoms
      integer, intent(in) :: n_freqs                 !! Number of frequencies
      type(thermochemistry_result_t), intent(out) :: result  !! Output results
      real(dp), intent(in), optional :: temperature  !! Temperature in K (default 298.15)
      real(dp), intent(in), optional :: pressure     !! Pressure in atm (default 1.0)
      integer, intent(in), optional :: symmetry_number    !! Symmetry number (default 1)
      integer, intent(in), optional :: spin_multiplicity  !! Spin multiplicity (default 1)

      real(dp) :: center_of_mass(3)
      real(dp) :: principal_axes(3, 3)
      real(dp) :: T, P
      integer :: sigma, mult
      real(dp) :: S_total

      ! Set parameters
      T = DEFAULT_TEMPERATURE
      P = DEFAULT_PRESSURE
      sigma = DEFAULT_SYMMETRY_NUMBER
      mult = DEFAULT_SPIN_MULTIPLICITY

      if (present(temperature)) T = temperature
      if (present(pressure)) P = pressure
      if (present(symmetry_number)) sigma = symmetry_number
      if (present(spin_multiplicity)) mult = spin_multiplicity

      result%temperature = T
      result%pressure = P
      result%symmetry_number = sigma
      result%spin_multiplicity = mult

      ! Compute moments of inertia
      call compute_moments_of_inertia(coords, atomic_numbers, n_atoms, &
                                      center_of_mass, result%moments, principal_axes, &
                                      result%is_linear, result%total_mass)

      ! Compute rotational constants
      call compute_rotational_constants(result%moments, result%is_linear, result%rot_const)

      ! Compute ZPE
      call compute_zpe(frequencies, n_freqs, result%n_real_freqs, &
                       result%zpe_hartree, result%zpe_kcalmol)
      result%n_imag_freqs = n_freqs - result%n_real_freqs

      ! Compute translational contributions
      call compute_translational_thermo(result%total_mass, T, P, &
                                        result%E_trans, result%S_trans, result%Cv_trans)

      ! Compute rotational contributions
      call compute_rotational_thermo(result%moments, T, sigma, result%is_linear, &
                                     result%E_rot, result%S_rot, result%Cv_rot)

      ! Compute vibrational contributions (thermal, excluding ZPE)
      call compute_vibrational_thermo(frequencies, n_freqs, T, &
                                      result%E_vib, result%S_vib, result%Cv_vib)

      ! Compute electronic contributions
      result%E_elec = 0.0_dp  ! Ground state only
      call compute_electronic_entropy(mult, result%S_elec)

      ! Compute partition functions
      call compute_partition_functions(result%total_mass, result%moments, frequencies, n_freqs, &
                                       T, P, sigma, result%is_linear, &
                                       result%q_trans, result%q_rot, result%q_vib)

      ! Compute summary quantities
      ! Thermal correction to energy = ZPE + E_trans + E_rot + E_vib + E_elec
      result%thermal_correction_energy = result%zpe_hartree + &
                                         result%E_trans + result%E_rot + result%E_vib + result%E_elec

      ! Thermal correction to enthalpy = E_tot + RT
      result%thermal_correction_enthalpy = result%thermal_correction_energy + R_HARTREE*T

      ! Total entropy in Hartree/K
      S_total = (result%S_trans + result%S_rot + result%S_vib + result%S_elec)/HARTREE_TO_CALMOL

      ! Thermal correction to Gibbs = H - TS
      result%thermal_correction_gibbs = result%thermal_correction_enthalpy - T*S_total

   end subroutine compute_thermochemistry