1!--------------------------------------------------------------------------------------------------!
2!   CP2K: A general program to perform molecular dynamics simulations                              !
3!   Copyright (C) 2000 - 2019  CP2K developers group                                               !
4!--------------------------------------------------------------------------------------------------!
5
6! **************************************************************************************************
7!> \brief Definition of physical constants:
8!>
9!>      a_bohr      : Bohr radius [m]
10!>      a_fine      : Fine-structure constant
11!>      a_mass      : Atomic mass unit [kg]; conversion factor [u] -> [kg]
12!>      angstrom    : Conversion factor [Bohr] -> [Angstrom]
13!>      bar         : Conversion factor [a.u.] -> [bar]
14!>      bohr        : Conversion factor [Angstrom] -> [Bohr]
15!>      boltzmann   : Boltzmann constant [J/K]
16!>      c_light     : Speed of light in vacuum [m/s]
17!>      c_light_au  : Speed of light in vacuum [a.u.]
18!>      e_charge    : Elementary charge [C]
19!>      e_mass      : Electron mass [kg]
20!>      e_gfactor   : Electron g factor [ ]
21!>      esu         : Conversion factors [a.u.] -> [esu]
22!>      evolt       : Conversion factor [a.u.] -> [eV]
23!>      femtoseconds: Conversion factor [a.u.] -> [fs]
24!>      h_bar       : Planck constant [J*s]
25!>      h_planck    : Planck constant [J*s]
26!>      hertz       : Conversion factor [a.u.] -> [Hz]
27!>      joule       : Conversion factor [a.u.] -> [J]
28!>      kcalmol     : Conversion factor [a.u.] -> [kcal/mol]
29!>      kelvin      : Conversion factor [a.u.] -> [K]
30!>      kjmol       : Conversion factor [a.u.] -> [kJ/mol]
31!>      massunit    : Conversion factor [u] -> [a.u.]
32!>      mu_perm     : Magnetic constant or permeability of vacuum [N/A**2]
33!>      n_avogadro  : Avogadro constant [1/mol]
34!>      newton      : Conversion factor [a.u.] -> [N]
35!>      pascal      : Conversion factor [a.u.] -> [Pa]
36!>      permittivity: Electric constant or permittivity of vacuum [F/m]
37!>      picoseconds : Conversion factor [a.u.] -> [ps]
38!>      rydberg     : Rydberg constant [1/m]
39!>      seconds     : Conversion factor [a.u.] -> [s]
40!>      vibfac      : Conversion factor [a.u./Bohr**2] -> [1/cm]
41!>      wavenumbers : Conversion factor [a.u.] -> [1/cm]
42!>      debye       : Conversion factor [a.u.] -> Debye
43!> \note
44!>      Fundamental physical constants (SI units)
45!>      Literature: - P. J. Mohr and B. N. Taylor,
46!>                    "CODATA recommended values of the fundamental physical
47!>                     constants: 1998 Rev. Mod. Phys. 72, 351-495 (2000)
48!>                  - P. J. Mohr and B. N. Taylor,
49!>                    "CODATA recommended values of the fundamental physical
50!>                     constants: 2002", Rev. Mod. Phys. 77, 1 (2005).
51!>                  - P. J. Mohr, B. N. Taylor, and D. B. Newell,
52!>                    "CODATA recommended values of the fundamental physical
53!>                     constants: 2006 Rev. Mod. Phys. 80, 633 (2008)
54!>                  - P. J. Mohr, B. N. Taylor, and D. B. Newell,
55!>                    "CODATA recommended values of the fundamental physical
56!>                     constants: 2010", Rev. Mod. Phys. 84, 1527-1605 (2012)
57!> \par History
58!>      - Adapted for use in CP2K (JGH)
59!>      - Updated to CODATA 1998 and cleaned (05.09.2003,MK)
60!>      - Updated to CODATA 2006. (26.03.2008,AK)
61!>      - Updated to CODATA 2010. (10.12.2012,MK)
62!>      - Turned constants into Fortran parameters (2014, Ole Schuett)
63!>      - Remove all but CODATA 2006 (2015, Ole Schuett)
64!> \author Matthias Krack (MK)
65! **************************************************************************************************
66MODULE physcon
67
68   USE kinds,                           ONLY: dp
69   USE mathconstants,                   ONLY: pi
70
71   IMPLICIT NONE
72
73   PRIVATE
74
75   PUBLIC :: a_bohr, a_fine, a_mass, angstrom, atm, bar, bohr, boltzmann, c_light_au, &
76             debye, e_charge, e_gfactor, e_mass, evolt, femtoseconds, h_bar, &
77             hertz, joule, kcalmol, kelvin, kjmol, massunit, mu_perm, n_avogadro, newton, &
78             p_mass, pascal, picoseconds, seconds, vibfac, wavenumbers
79
80   PUBLIC :: write_physcon
81
82   ! CP2K uses the CODATA from 2006
83
84   ! Exact constants
85   ! Speed of light in vacuum [m/s]
86   REAL(KIND=dp), PARAMETER :: c_light = 299792458.0_dp
87   ! Speed of light in vacuum, in atomic units (=1/a_fine)
88   REAL(KIND=dp), PARAMETER :: c_light_au = 137.035999679_dp
89
90   ! Magnetic constant or permeability of vacuum [N/A**2]
91   REAL(KIND=dp), PARAMETER :: mu_perm = 4.0_dp*pi*1.0E-7_dp
92
93   ! Electric constant or permittivity of vacuum [F/m]
94   REAL(KIND=dp), PARAMETER :: permittivity = 1.0_dp/(mu_perm*c_light**2)
95
96   ! Recommended fundamental constants of physics
97   ! and chemistry based on the 2006 adjustment
98
99   ! Planck constant [J*s]
100   REAL(KIND=dp), PARAMETER :: h_planck = 6.62606896E-34_dp
101   REAL(KIND=dp), PARAMETER :: h_bar = h_planck/(2.0_dp*pi)
102
103   ! Elementary charge [C]
104   REAL(KIND=dp), PARAMETER :: e_charge = 1.602176487E-19_dp
105
106   ! Electron mass [kg]
107   REAL(KIND=dp), PARAMETER :: e_mass = 9.10938215E-31_dp
108
109   ! Proton mass [kg]
110   REAL(KIND=dp), PARAMETER :: p_mass = 1.672621637E-27_dp
111
112   ! Electron g factor [ ]
113   REAL(KIND=dp), PARAMETER :: e_gfactor = -2.0023193043622_dp
114
115   ! Fine-structure constant
116!MK a_fine = 0.5_dp*mu_perm*c_light*e_charge**2/h_planck
117   REAL(KIND=dp), PARAMETER :: a_fine = 7.2973525376E-3_dp
118
119   ! Rydberg constant [1/m]
120!MK rydberg = 0.5_dp*e_mass*c_light*a_fine**2/h_planck
121   REAL(KIND=dp), PARAMETER :: rydberg = 10973731.568527_dp
122
123   ! Avogadro constant [1/mol]
124   REAL(KIND=dp), PARAMETER :: n_avogadro = 6.02214179E+23_dp
125
126   ! Boltzmann constant [J/K]
127   REAL(KIND=dp), PARAMETER :: boltzmann = 1.3806504E-23_dp
128
129   ! Atomic mass unit [kg]; conversion factor [u] -> [kg]
130   REAL(KIND=dp), PARAMETER :: a_mass = 1.660538782E-27_dp
131
132   ! Bohr radius [m]
133!MK a_bohr = a_fine/(4.0_dp*pi*rydberg)
134   REAL(KIND=dp), PARAMETER :: a_bohr = 0.52917720859E-10_dp
135
136   ! Conversion factors
137
138   ! [u] -> [a.u.]
139   REAL(KIND=dp), PARAMETER :: massunit = a_mass/e_mass
140
141   ! [Bohr] -> [Angstrom]
142   REAL(KIND=dp), PARAMETER :: angstrom = 1.0E+10_dp*a_bohr
143
144   ! [Angstrom] -> [Bohr]
145   REAL(KIND=dp), PARAMETER :: bohr = 1.0_dp/angstrom
146
147   ! [a.u.] -> [s]
148   REAL(KIND=dp), PARAMETER :: seconds = 1.0_dp/(4.0_dp*pi*rydberg*c_light)
149
150   ! [a.u.] -> [fs]
151   REAL(KIND=dp), PARAMETER :: femtoseconds = 1.0E+15_dp*seconds
152
153   ! [a.u.] -> [ps]
154   REAL(KIND=dp), PARAMETER :: picoseconds = 1.0E+12_dp*seconds
155
156   ! [a.u.] -> [J]
157   REAL(KIND=dp), PARAMETER :: joule = 2.0_dp*rydberg*h_planck*c_light
158
159   ! [a.u.] -> [N]
160   REAL(KIND=dp), PARAMETER :: newton = joule/a_bohr
161
162   ! [a.u.] -> [K]
163   REAL(KIND=dp), PARAMETER :: kelvin = joule/boltzmann
164
165   ! [a.u.] -> [kJ/mol]
166   REAL(KIND=dp), PARAMETER :: kjmol = 0.001_dp*joule*n_avogadro
167
168   ! [a.u.] -> [kcal/mol]
169   REAL(KIND=dp), PARAMETER :: kcalmol = kjmol/4.184_dp
170
171   ! [a.u.] -> [Pa]
172   REAL(KIND=dp), PARAMETER :: pascal = joule/a_bohr**3
173
174   ! [a.u.] -> [bar]
175   REAL(KIND=dp), PARAMETER :: bar = pascal/1.0E+5_dp
176
177   ! [a.u.] -> [atm]
178   REAL(KIND=dp), PARAMETER :: atm = pascal/1.013250E+5_dp
179
180   ! [a.u.] -> [eV]
181   REAL(KIND=dp), PARAMETER :: evolt = joule/e_charge
182
183   ! [a.u.] -> [Hz]
184   REAL(KIND=dp), PARAMETER :: hertz = joule/h_planck
185
186   ! [a.u./Bohr**2] -> [1/cm] (wave numbers)
187   REAL(KIND=dp), PARAMETER :: vibfac = 5.0_dp*SQRT(kjmol)/(pi*a_bohr*c_light)
188
189   ! [a.u.] -> [1/cm] (wave numbers)
190   REAL(KIND=dp), PARAMETER :: wavenumbers = 0.02_dp*rydberg
191
192   ! [a.u.] -> [esu] (electrostatic units)
193   REAL(KIND=dp), PARAMETER :: esu_1 = 1.0E+21_dp*a_bohr*c_light*e_charge
194   !REAL(KIND=dp), PARAMETER :: esu_2  = esu_1/bohr
195   !REAL(KIND=dp), PARAMETER :: esu_3  = esu_2/bohr
196   !REAL(KIND=dp), PARAMETER :: esu(3) = (/ esu_1, esu_2, esu_3 /)
197
198   ! [a.u.] -> [debye] (electrostatic units)
199   REAL(KIND=dp), PARAMETER :: debye = esu_1
200
201CONTAINS
202
203! **************************************************************************************************
204!> \brief  Write all basic physical constants used by CP2K to a logical
205!>           output unit.
206!> \param output_unit ...
207!> \date    14.11.2000
208!> \par History
209!>       - Updated to CODATA 1998 and cleaned (05.09.2003,MK)
210!>       - Updated to CODATA 2006. (26.03.2008,AK)
211!> \author  JGH
212!> \version 1.1
213! **************************************************************************************************
214   SUBROUTINE write_physcon(output_unit)
215
216      INTEGER, INTENT(IN)                                :: output_unit
217
218      WRITE (UNIT=output_unit, FMT="(T2,/,T2,A,/,/,(T2,A))") &
219         "*** Fundamental physical constants (SI units) ***", &
220         "*** Literature: B. J. Mohr and B. N. Taylor,", &
221         "***             CODATA recommended values of the fundamental physical", &
222         "***             constants: 2006, Web Version 5.1", &
223         "***             http://physics.nist.gov/constants"
224
225      WRITE (UNIT=output_unit, FMT="(/,T2,A,T61,ES20.14)") &
226         "Speed of light in vacuum [m/s]", c_light
227      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
228         "Magnetic constant or permeability of vacuum [N/A**2]", mu_perm
229      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
230         "Electric constant or permittivity of vacuum [F/m]", permittivity
231      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
232         "Planck constant (h) [J*s]", h_planck
233      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
234         "Planck constant (h-bar) [J*s]", h_bar
235      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
236         "Elementary charge [C]", e_charge
237      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
238         "Electron mass [kg]", e_mass
239      WRITE (UNIT=output_unit, FMT="(T2,A,T60,ES21.14)") &
240         "Electron g factor [ ]", e_gfactor
241      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
242         "Proton mass [kg]", p_mass
243      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
244         "Fine-structure constant", a_fine
245      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
246         "Rydberg constant [1/m]", rydberg
247      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
248         "Avogadro constant [1/mol]", n_avogadro
249      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
250         "Boltzmann constant [J/K]", boltzmann
251      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
252         "Atomic mass unit [kg]", a_mass
253      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
254         "Bohr radius [m]", a_bohr
255
256      ! Conversion factors
257
258      WRITE (UNIT=output_unit, FMT="(/,T2,A,/)") &
259         "*** Conversion factors ***"
260
261      WRITE (UNIT=output_unit, FMT="(T2,A,T61,ES20.14)") &
262         "[u] -> [a.u.]", massunit, &
263         "[Angstrom] -> [Bohr] = [a.u.]", bohr, &
264         "[a.u.] = [Bohr] -> [Angstrom]", angstrom, &
265         "[a.u.] -> [s]", seconds, &
266         "[a.u.] -> [fs]", femtoseconds, &
267         "[a.u.] -> [J]", joule, &
268         "[a.u.] -> [N]", newton, &
269         "[a.u.] -> [K]", kelvin, &
270         "[a.u.] -> [kJ/mol]", kjmol, &
271         "[a.u.] -> [kcal/mol]", kcalmol, &
272         "[a.u.] -> [Pa]", pascal, &
273         "[a.u.] -> [bar]", bar, &
274         "[a.u.] -> [atm]", atm, &
275         "[a.u.] -> [eV]", evolt, &
276         "[a.u.] -> [Hz]", hertz, &
277         "[a.u.] -> [1/cm] (wave numbers)", wavenumbers, &
278         "[a.u./Bohr**2] -> [1/cm]", vibfac
279      WRITE (UNIT=output_unit, FMT="(T2,A)") ""
280
281   END SUBROUTINE write_physcon
282
283END MODULE physcon
284