ARTS  2.4.0(git:4fb77825)
constants.h
Go to the documentation of this file.
1 /* Copyright (C) 2019
2  * Richard Larsson <ric.larsson@gmail.com>
3  *
4  * This program is free software; you can redistribute it and/or modify it
5  * under the terms of the GNU General Public License as published by the
6  * Free Software Foundation; either version 2, or (at your option) any
7  * later version.
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program; if not, write to the Free Software
16  * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
17  * USA. */
18 
54 #ifndef CONSTANTS_IN_ARTS_H
55 #define CONSTANTS_IN_ARTS_H
56 
57 #include <cmath>
58 #include "matpack.h"
59 
61 namespace Constant {
63 template <class T>
64 constexpr T pow2(T x) {
65  return x * x;
66 }
67 
69 template <class T>
70 constexpr T pow3(T x) {
71  return pow2(x) * x;
72 }
73 
75 template <class T>
76 constexpr T pow4(T x) {
77  return pow2(pow2(x));
78 }
79 
111 static constexpr Numeric pi =
112  3.141592653589793238462643383279502884197169399375105820974944592307816406286;
113 
115 static constexpr Numeric inv_pi =
116  0.3183098861837906715377675267450287240689192914809128974953346881177935952685;
117 
119 static constexpr Numeric two_pi =
120  6.283185307179586476925286766559005768394338798750211641949889184615632812572;
121 
123 static constexpr Numeric inv_two_pi =
124  0.1591549430918953357688837633725143620344596457404564487476673440588967976342;
125 
127 static constexpr Numeric sqrt_pi =
128  1.772453850905516027298167483341145182797549456122387128213807789852911284591;
129 
131 static constexpr Numeric inv_sqrt_pi =
132  0.5641895835477562869480794515607725858440506293289988568440857217106424684415;
133 
135 static constexpr Numeric euler =
136  2.718281828459045235360287471352662497757247093699959574966967627724076630354;
137 
139 static constexpr Numeric inv_euler =
140  0.3678794411714423215955237701614608674458111310317678345078368016974614957448;
141 
143 static constexpr Numeric log10_euler =
144  0.4342944819032518276511289189166050822943970058036665661144537831658646492089;
145 
147 static constexpr Numeric ln_10 =
148  2.302585092994045684017991454684364207601101488628772976033327900967572609677;
149 
151 static constexpr Numeric sqrt_2 =
152  1.414213562373095048801688724209698078569671875376948073176679737990732478462;
153 
155 static constexpr Numeric inv_sqrt_2 =
156  0.7071067811865475244008443621048490392848359376884740365883398689953662392311;
157 
159 static constexpr Numeric ln_2 =
160  0.6931471805599453094172321214581765680755001343602552541206800094933936219697;
161 
163 static constexpr Numeric inv_ln_2 =
164  1.442695040888963407359924681001892137426645954152985934135449406931109219181;
165 
167 static constexpr Numeric sqrt_ln_2 =
168  0.8325546111576977563531646448952010476305888522644407291668291172340794351973;
169 
171 static constexpr Numeric inv_sqrt_ln_2 =
172  1.201122408786449794857803286095221722566764028068699423868879896733837175546;
173 
178 static constexpr Numeric Delta_nu_Cs = 9192631770;
179 
184 static constexpr Numeric speed_of_light = 299792458;
185 
187 static constexpr Numeric c = speed_of_light;
188 
193 static constexpr Numeric planck_constant = 6.62607015e-34;
194 
196 static constexpr Numeric h = planck_constant;
197 
199 static constexpr Numeric reduced_planck_constant = h * inv_two_pi;
200 
202 static constexpr Numeric h_bar = reduced_planck_constant;
203 
208 static constexpr Numeric elementary_charge = 1.602176634e-19;
209 
211 static constexpr Numeric e = elementary_charge;
212 
217 static constexpr Numeric boltzmann_constant = 1.380649e-23;
218 
220 static constexpr Numeric k = boltzmann_constant;
221 
226 static constexpr Numeric avogadro_constant = 6.02214076e23;
227 
229 static constexpr Numeric NA = avogadro_constant;
230 
235 static constexpr Numeric K_cd = 683;
236 
242 static constexpr Numeric fine_structure_constant = 7.2973525693e-3;
243 
245 static constexpr Numeric alpha = fine_structure_constant;
246 
252 static constexpr Numeric rydberg_constant = 10973731.568160;
253 
255 static constexpr Numeric R_inf = rydberg_constant;
256 
258 static constexpr Numeric magnetic_constant = 2 * h * alpha / (c * pow2(e));
259 
261 static constexpr Numeric mu_0 = magnetic_constant;
262 
264 static constexpr Numeric vacuum_permittivity = pow2(e) / (2 * h * c * alpha);
265 
267 static constexpr Numeric epsilon_0 = vacuum_permittivity;
268 
270 static constexpr Numeric electron_mass = 2 * h * R_inf / (c * pow2(alpha));
271 
273 static constexpr Numeric m_e = electron_mass;
274 
280 static constexpr Numeric unified_atomic_mass_unit = 1.66053906660e-27;
281 
283 static constexpr Numeric m_u = unified_atomic_mass_unit;
284 
290 static constexpr Numeric mass_ratio_electrons_per_proton = 1'836.152'673'43;
291 
293 static constexpr Numeric proton_mass = electron_mass * mass_ratio_electrons_per_proton;
294 
300 static constexpr Numeric mass_ratio_electrons_per_neutron = 1'838.683'661'73;
301 
303 static constexpr Numeric neutron_mass = electron_mass * mass_ratio_electrons_per_neutron;
304 
306 static constexpr Numeric bohr_magneton = e * h_bar / (2 * m_e);
307 
309 static constexpr Numeric ideal_gas_constant = k * NA;
310 
312 static constexpr Numeric R = ideal_gas_constant;
313 
315 static constexpr Numeric doppler_broadening_const_squared = 2000 * R / pow2(c);
316 }; // namespace Constant
317 
319 namespace Conversion {
320 using namespace Constant;
321 
323 static constexpr Numeric DEG2RAD = pi / 180;
324 static constexpr Numeric RAD2DEG = 1 / DEG2RAD;
326 template <class T>
327 constexpr Numeric deg2rad(T x) {
328  return x * DEG2RAD;
329 }
330 
332 template <class T>
333 constexpr Numeric rad2deg(T x) {
334  return x * RAD2DEG;
335 }
336 
338 template <class T>
340  return std::cos(deg2rad(x));
341 }
342 
344 template <class T>
346  return std::sin(deg2rad(x));
347 }
348 
350 template <class T>
352  return std::tan(deg2rad(x));
353 }
354 
356 template <class T>
358  return rad2deg(std::acos(x));
359 }
360 
362 template <class T>
364  return rad2deg(std::asin(x));
365 }
366 
368 template <class T>
370  return rad2deg(std::atan(x));
371 }
372 
374 template <class T1, class T2>
375 Numeric atan2d(T1 y, T2 x) {
376  return rad2deg(std::atan2(y, x));
377 }
378 
380 static constexpr Numeric KAYCM2FREQ = 100 * c;
381 static constexpr Numeric FREQ2KAYCM = 1 / KAYCM2FREQ;
382 template <class T>
383 constexpr Numeric kaycm2freq(T x) {
384  return x * KAYCM2FREQ;
385 }
386 template <class T>
387 constexpr Numeric freq2kaycm(T x) {
388  return x * FREQ2KAYCM;
389 }
390 
392 static constexpr Numeric ANGCM2FREQ = KAYCM2FREQ * inv_two_pi;
393 static constexpr Numeric FREQ2ANGCM = 1 / ANGCM2FREQ;
394 template <class T>
395 constexpr Numeric angcm2freq(T x) {
396  return x * ANGCM2FREQ;
397 }
398 template <class T>
399 constexpr Numeric freq2angcm(T x) {
400  return x * FREQ2ANGCM;
401 }
402 
404 template <class T>
405 constexpr Numeric angfreq2freq(T x) {
406  return x * inv_two_pi;
407 }
408 template <class T>
409 constexpr Numeric freq2angfreq(T x) {
410  return x * two_pi;
411 }
412 
414 template <class T>
415 constexpr Numeric wavelen2freq(T x) {
416  return c / x;
417 }
418 template <class T>
419 constexpr Numeric freq2wavelen(T x) {
420  return c / x;
421 }
422 
424 static constexpr Numeric ATM2PA = 101325;
425 static constexpr Numeric PA2ATM = 1 / ATM2PA;
426 template <class T>
427 constexpr Numeric atm2pa(T x) {
428  return x * ATM2PA;
429 }
430 template <class T>
431 constexpr Numeric pa2atm(T x) {
432  return x * PA2ATM;
433 }
434 
436 static constexpr Numeric TORR2PA = ATM2PA / 760;
437 static constexpr Numeric PA2TORR = 1 / TORR2PA;
438 template <class T>
439 constexpr Numeric torr2pa(T x) {
440  return x * TORR2PA;
441 }
442 template <class T>
443 constexpr Numeric pa2torr(T x) {
444  return x * PA2TORR;
445 }
446 
448 static constexpr Numeric CEL2KEL = 273.15;
449 template <class T>
450 constexpr Numeric celsius2kelvin(T x) {
451  return x + CEL2KEL;
452 }
453 template <class T>
454 constexpr Numeric kelvin2celsius(T x) {
455  return x - CEL2KEL;
456 }
457 
459 static constexpr Numeric HITRAN2ARTS_LS = KAYCM2FREQ * 1e-4;
460 static constexpr Numeric ARTS2HITRAN_LS = 1 / HITRAN2ARTS_LS;
461 template <class T>
462 constexpr T hitran2arts_linestrength(T x) {
463  return x * HITRAN2ARTS_LS;
464 }
465 template <class T>
466 constexpr T arts2hitran_linestrength(T x) {
467  return x * ARTS2HITRAN_LS;
468 }
469 
471 static constexpr Numeric HITRAN2ARTS_GAMMA = KAYCM2FREQ / ATM2PA;
472 static constexpr Numeric ARTS2HITRAN_GAMMA = 1 / HITRAN2ARTS_GAMMA;
473 template <class T>
474 constexpr T hitran2arts_broadening(T x) {
475  return x * HITRAN2ARTS_GAMMA;
476 }
477 template <class T>
478 constexpr T arts2hitran_broadening(T x) {
479  return x * ARTS2HITRAN_GAMMA;
480 }
481 
483 static constexpr Numeric HITRAN2ARTS_ENERGY = h * KAYCM2FREQ;
484 static constexpr Numeric ARTS2HITRAN_ENERGY = 1 / HITRAN2ARTS_ENERGY;
485 template <class T>
486 constexpr T hitran2arts_energy(T x) {
487  return x * HITRAN2ARTS_ENERGY;
488 }
489 template <class T>
490 constexpr T arts2hitran_energy(T x) {
491  return x * ARTS2HITRAN_ENERGY;
492 }
493 
494 template <class T>
495 constexpr T angstrom2meter(T x) {
496  return x * 1e-10;
497 }
498 template <class T>
499 constexpr T meter2angstrom(T x) {
500  return x * 1e10;
501 }
502 }; // namespace Conversion
503 
504 #endif
Conversion::freq2angcm
constexpr Numeric freq2angcm(T x)
Definition: constants.h:399
matpack.h
Conversion::kelvin2celsius
constexpr Numeric kelvin2celsius(T x)
Definition: constants.h:454
Conversion::tand
Numeric tand(T x)
Returns the tangent of the deg2rad of the input
Definition: constants.h:351
Conversion::torr2pa
constexpr Numeric torr2pa(T x)
Definition: constants.h:439
Conversion::arts2hitran_broadening
constexpr T arts2hitran_broadening(T x)
Definition: constants.h:478
Conversion::kaycm2freq
constexpr Numeric kaycm2freq(T x)
Definition: constants.h:383
Conversion::atand
Numeric atand(T x)
Returns rad2deg of the arc-tangent of the input
Definition: constants.h:369
Conversion::freq2angfreq
constexpr Numeric freq2angfreq(T x)
Definition: constants.h:409
RAD2DEG
const Numeric RAD2DEG
Global constant, conversion from radians to degrees.
Conversion::acosd
Numeric acosd(T x)
Returns rad2deg of the arc-cosine of the input
Definition: constants.h:357
Conversion::freq2kaycm
constexpr Numeric freq2kaycm(T x)
Definition: constants.h:387
ARTS::Var::y
Vector y(Workspace &ws) noexcept
Definition: autoarts.h:7401
DEG2RAD
const Numeric DEG2RAD
Global constant, conversion from degrees to radians.
Conversion::asind
Numeric asind(T x)
Returns rad2deg of the arc-sine of the input
Definition: constants.h:363
Conversion::arts2hitran_energy
constexpr T arts2hitran_energy(T x)
Definition: constants.h:490
Conversion::wavelen2freq
constexpr Numeric wavelen2freq(T x)
Conversion wavelength to frequency and back.
Definition: constants.h:415
Conversion::angstrom2meter
constexpr T angstrom2meter(T x)
Definition: constants.h:495
Constant::pow4
constexpr T pow4(T x)
power of four
Definition: constants.h:76
Conversion::hitran2arts_linestrength
constexpr T hitran2arts_linestrength(T x)
Definition: constants.h:462
TORR2PA
const Numeric TORR2PA
Global constant, converts torr to Pa.
Conversion::celsius2kelvin
constexpr Numeric celsius2kelvin(T x)
Definition: constants.h:450
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Conversion::freq2wavelen
constexpr Numeric freq2wavelen(T x)
Definition: constants.h:419
Conversion::meter2angstrom
constexpr T meter2angstrom(T x)
Definition: constants.h:499
ATM2PA
const Numeric ATM2PA
Global constant, converts atm to Pa.
Conversion::pa2torr
constexpr Numeric pa2torr(T x)
Definition: constants.h:443
Conversion::hitran2arts_energy
constexpr T hitran2arts_energy(T x)
Definition: constants.h:486
Conversion::hitran2arts_broadening
constexpr T hitran2arts_broadening(T x)
Definition: constants.h:474
Conversion::pa2atm
constexpr Numeric pa2atm(T x)
Definition: constants.h:431
Conversion::rad2deg
constexpr Numeric rad2deg(T x)
Converts radians to degrees
Definition: constants.h:333
Constant::pow3
constexpr T pow3(T x)
power of three
Definition: constants.h:70
Conversion::arts2hitran_linestrength
constexpr T arts2hitran_linestrength(T x)
Definition: constants.h:466
Conversion::deg2rad
constexpr Numeric deg2rad(T x)
Converts degrees to radians
Definition: constants.h:327
Conversion::atm2pa
constexpr Numeric atm2pa(T x)
Definition: constants.h:427
Conversion::cosd
Numeric cosd(T x)
Returns the cosine of the deg2rad of the input
Definition: constants.h:339
Constant::pow2
constexpr T pow2(T x)
power of two
Definition: constants.h:64
ARTS::Var::x
Vector x(Workspace &ws) noexcept
Definition: autoarts.h:7346
Conversion
Namespace containing several practical unit conversions, physical and mathematical.
Definition: constants.h:319
Conversion::atan2d
Numeric atan2d(T1 y, T2 x)
Returns rad2deg of the arc-tangent of inputs #T1/#T2
Definition: constants.h:375
Conversion::sind
Numeric sind(T x)
Returns the sine of the deg2rad of the input
Definition: constants.h:345
Conversion::angfreq2freq
constexpr Numeric angfreq2freq(T x)
Conversion constant Angular frequency to frequency and back.
Definition: constants.h:405
Conversion::angcm2freq
constexpr Numeric angcm2freq(T x)
Definition: constants.h:395
Constant
Namespace containing several constants, physical and mathematical.
Definition: constants.h:61