ARTS 2.5.0 (git: 9ee3ac6c)
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 "enums.h"
58#include "matpack.h"
59#include <cmath>
60
62namespace Constant {
64template <class T>
65constexpr auto pow2(T x) -> decltype(x * x) {
66 return x * x;
67}
68
70template <class T>
71constexpr auto pow3(T x) -> decltype(pow2(x) * x) {
72 return pow2(x) * x;
73}
74
76template <class T>
77constexpr auto pow4(T x) -> decltype(pow2(pow2(x))) {
78 return pow2(pow2(x));
79}
80
112static constexpr Numeric pi =
113 3.141592653589793238462643383279502884197169399375105820974944592307816406286;
114
116static constexpr Numeric inv_pi =
117 0.3183098861837906715377675267450287240689192914809128974953346881177935952685;
118
120static constexpr Numeric two_pi =
121 6.283185307179586476925286766559005768394338798750211641949889184615632812572;
122
124static constexpr Numeric inv_two_pi =
125 0.1591549430918953357688837633725143620344596457404564487476673440588967976342;
126
128static constexpr Numeric sqrt_pi =
129 1.772453850905516027298167483341145182797549456122387128213807789852911284591;
130
132static constexpr Numeric inv_sqrt_pi =
133 0.5641895835477562869480794515607725858440506293289988568440857217106424684415;
134
136static constexpr Numeric euler =
137 2.718281828459045235360287471352662497757247093699959574966967627724076630354;
138
140static constexpr Numeric inv_euler =
141 0.3678794411714423215955237701614608674458111310317678345078368016974614957448;
142
144static constexpr Numeric log10_euler =
145 0.4342944819032518276511289189166050822943970058036665661144537831658646492089;
146
148static constexpr Numeric ln_10 =
149 2.302585092994045684017991454684364207601101488628772976033327900967572609677;
150
152static constexpr Numeric sqrt_2 =
153 1.414213562373095048801688724209698078569671875376948073176679737990732478462;
154
156static constexpr Numeric inv_sqrt_2 =
157 0.7071067811865475244008443621048490392848359376884740365883398689953662392311;
158
160static constexpr Numeric ln_2 =
161 0.6931471805599453094172321214581765680755001343602552541206800094933936219697;
162
164static constexpr Numeric inv_ln_2 =
165 1.442695040888963407359924681001892137426645954152985934135449406931109219181;
166
168static constexpr Numeric sqrt_ln_2 =
169 0.8325546111576977563531646448952010476305888522644407291668291172340794351973;
170
172static constexpr Numeric inv_sqrt_ln_2 =
173 1.201122408786449794857803286095221722566764028068699423868879896733837175546;
174
179static constexpr Numeric Delta_nu_Cs = 9192631770;
180
185static constexpr Numeric speed_of_light = 299792458;
186
188static constexpr Numeric c = speed_of_light;
189
194static constexpr Numeric planck_constant = 6.62607015e-34;
195
197static constexpr Numeric h = planck_constant;
198
200static constexpr Numeric reduced_planck_constant = h * inv_two_pi;
201
203static constexpr Numeric h_bar = reduced_planck_constant;
204
209static constexpr Numeric elementary_charge = 1.602176634e-19;
210
212static constexpr Numeric e = elementary_charge;
213
218static constexpr Numeric boltzmann_constant = 1.380649e-23;
219
221static constexpr Numeric k = boltzmann_constant;
222
227static constexpr Numeric avogadro_constant = 6.02214076e23;
228
230static constexpr Numeric NA = avogadro_constant;
231
236static constexpr Numeric K_cd = 683;
237
243static constexpr Numeric fine_structure_constant = 7.2973525693e-3;
244
246static constexpr Numeric alpha = fine_structure_constant;
247
253static constexpr Numeric rydberg_constant = 10973731.568160;
254
256static constexpr Numeric R_inf = rydberg_constant;
257
259static constexpr Numeric magnetic_constant = 2 * h * alpha / (c * pow2(e));
260
262static constexpr Numeric mu_0 = magnetic_constant;
263
265static constexpr Numeric vacuum_permittivity = pow2(e) / (2 * h * c * alpha);
266
268static constexpr Numeric epsilon_0 = vacuum_permittivity;
269
271static constexpr Numeric electron_mass = 2 * h * R_inf / (c * pow2(alpha));
272
274static constexpr Numeric m_e = electron_mass;
275
281static constexpr Numeric unified_atomic_mass_unit = 1.66053906660e-27;
282
284static constexpr Numeric m_u = unified_atomic_mass_unit;
285
291static constexpr Numeric mass_ratio_electrons_per_proton = 1'836.152'673'43;
292
294static constexpr Numeric proton_mass = electron_mass * mass_ratio_electrons_per_proton;
295
301static constexpr Numeric mass_ratio_electrons_per_neutron = 1'838.683'661'73;
302
304static constexpr Numeric neutron_mass = electron_mass * mass_ratio_electrons_per_neutron;
305
307static constexpr Numeric bohr_magneton = e * h_bar / (2 * m_e);
308
310static constexpr Numeric ideal_gas_constant = k * NA;
311
313static constexpr Numeric R = ideal_gas_constant;
314
316static constexpr Numeric doppler_broadening_const_squared = 2'000 * R / pow2(c);
317
319static constexpr Numeric one_degree_in_radians = pi / 180;
320}; // namespace Constant
321
323namespace Conversion {
324using namespace Constant;
325
327template <class T>
328constexpr auto deg2rad(T x) -> decltype(x * one_degree_in_radians) {
329 return x * one_degree_in_radians;
330}
331
333template <class T>
334constexpr auto rad2deg(T x) -> decltype(x / one_degree_in_radians) {
335 return x / one_degree_in_radians;
336}
337
339template <class T>
340auto cosd(T x) -> decltype(std::cos(deg2rad(x))) {
341 return std::cos(deg2rad(x));
342}
343
345template <class T>
346auto sind(T x) -> decltype(std::sin(deg2rad(x))) {
347 return std::sin(deg2rad(x));
348}
349
351template <class T>
352auto tand(T x) -> decltype(std::tan(deg2rad(x))) {
353 return std::tan(deg2rad(x));
354}
355
357template <class T>
358auto acosd(T x) -> decltype(rad2deg(std::acos(x))) {
359 return rad2deg(std::acos(x));
360}
361
363template <class T>
364auto asind(T x) -> decltype(rad2deg(std::asin(x))) {
365 return rad2deg(std::asin(x));
366}
367
369template <class T>
370auto atand(T x) -> decltype(rad2deg(std::atan(x))) {
371 return rad2deg(std::atan(x));
372}
373
375template <class T1, class T2>
376auto atan2d(T1 y, T2 x) -> decltype(rad2deg(std::atan2(y, x))) {
377 return rad2deg(std::atan2(y, x));
378}
379
381template <class T>
382constexpr auto kaycm2freq(T x) -> decltype(x * (100 * c)) {
383 return x * (100 * c);
384}
385
387template <class T>
388constexpr auto freq2kaycm(T x) -> decltype(x / (100 * c)) {
389 return x / (100 * c);
390}
391
393template <class T>
394constexpr auto angcm2freq(T x) -> decltype(kaycm2freq(inv_two_pi)) {
395 return x * kaycm2freq(inv_two_pi);
396}
397
399template <class T>
400constexpr auto freq2angcm(T x) -> decltype(x / kaycm2freq(inv_two_pi)) {
401 return x / kaycm2freq(inv_two_pi);
402}
403
405template <class T>
406constexpr auto angfreq2freq(T x) -> decltype(x * inv_two_pi) {
407 return x * inv_two_pi;
408}
409
411template <class T>
412constexpr auto freq2angfreq(T x) -> decltype(x * two_pi) {
413 return x * two_pi;
414}
415
417template <class T>
418constexpr auto wavelen2freq(T x) -> decltype(c / x) {
419 return c / x;
420}
421
423template <class T>
424constexpr auto freq2wavelen(T x) -> decltype(c / x) {
425 return c / x;
426}
427
429template <class T>
430constexpr auto atm2pa(T x) -> decltype(x * 101'325.0) {
431 return x * 101'325.0;
432}
433
435template <class T>
436constexpr auto pa2atm(T x) -> decltype(x / 101'325.0) {
437 return x / 101'325.0;
438}
439
441template <class T>
442constexpr auto torr2pa(T x) -> decltype(x * atm2pa(1.0 / 760.0)) {
443 return x * atm2pa(1.0 / 760.0);
444}
445
447template <class T>
448constexpr auto pa2torr(T x) -> decltype(x / atm2pa(1.0 / 760.0)) {
449 return x / atm2pa(1.0 / 760.0);
450}
451
453template <class T>
454constexpr auto mhz_per_torr2hz_per_pa(T x) -> decltype(x * pa2torr(1e6)) {
455 return x * pa2torr(1e6);
456}
457
459template <class T>
460constexpr auto celsius2kelvin(T x) -> decltype(x + 273.15) {
461 return x + 273.15;
462}
463
465template <class T>
466constexpr auto kelvin2celsius(T x) -> decltype(x - 273.15) {
467 return x - 273.15;
468}
469
471template <class T>
472constexpr auto kaycm_per_cmsquared2hz_per_msquared(T x) -> decltype(x * kaycm2freq(1e-4)) {
473 return x * kaycm2freq(1e-4);
474}
475
477template <class T>
478constexpr auto hz_per_msquared2kaycm_per_cmsquared(T x) -> decltype(x * freq2kaycm(1e4)) {
479 return x * freq2kaycm(1e4);
480}
481
483template <class T>
484constexpr auto kaycm_per_atm2hz_per_pa(T x) -> decltype(x * kaycm2freq(pa2atm(1))) {
485 return x * kaycm2freq(pa2atm(1));
486}
487
489template <class T>
490constexpr auto hz_per_pa2kaycm_per_atm(T x) -> decltype(x * freq2kaycm(atm2pa(1))) {
491 return x * freq2kaycm(atm2pa(1));
492}
493
495template <class T>
496constexpr auto kaycm2joule(T x) -> decltype(x * kaycm2freq(h)) {
497 return x * kaycm2freq(h);
498}
499
501template <class T>
502constexpr auto hz2joule(T x) -> decltype(x * h) {
503 return x * h;
504}
505
507template <class T>
508constexpr auto mhz2joule(T x) -> decltype(hz2joule(x) * 1e6) {
509 return hz2joule(x) * 1e6;
510}
511
513template <class T>
514constexpr auto kelvin2joule(T x) -> decltype(x * k) {
515 return x * k;
516}
517
519template <class T>
520constexpr auto joule2hz(T x) -> decltype(x / h) {
521 return x / h;
522}
523
525template <class T>
526constexpr auto joule2kaycm(T x) -> decltype(x / kaycm2freq(h)) {
527 return x / kaycm2freq(h);
528}
529
531template <class T>
532constexpr auto angstrom2meter(T x) -> decltype(x * 1e-10) {
533 return x * 1e-10;
534}
535
537template <class T>
538constexpr auto meter2angstrom(T x) -> decltype(x * 1e10) {
539 return x * 1e10;
540}
541}; // namespace Conversion
542
543namespace Options {
546 hour, hours, h,
547 minute, minutes, min,
548 second, seconds, s)
549
550
551 ENUMCLASS(LineShapeCoeff, char,
555 X3)
556
558 ENUMCLASS(WindMagJacobian, char,
559 u, v, w, strength)
560
562 ENUMCLASS(BasicCatParamJacobian, char,
564 LineCenter)
565
567 ENUMCLASS(HitranType, char,
568 Pre2004, // 2004 version changed the .par-length
569 From2004To2012, // New par length but old isotopologues order
570 Post2012, // 2012 version changed the order of isotopologues
571 Online // Onine expects a modern .par line followed by Upper then Lower quantum numbers
572 )
573
575 ENUMCLASS(LblSpeedup, char,
578 LinearIndependent)
579} // namespace Options
580
581#endif
std::chrono::duration< Numeric > TimeStep
A duration of time, 1 full tick should be 1 second.
Definition: artstime.h:38
#define min(a, b)
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Namespace containing several constants, physical and mathematical.
Definition: constants.h:62
constexpr auto pow3(T x) -> decltype(pow2(x) *x)
power of three
Definition: constants.h:71
constexpr auto pow4(T x) -> decltype(pow2(pow2(x)))
power of four
Definition: constants.h:77
constexpr auto pow2(T x) -> decltype(x *x)
power of two
Definition: constants.h:65
Namespace containing several practical unit conversions, physical and mathematical.
Definition: constants.h:323
constexpr auto freq2angfreq(T x) -> decltype(x *two_pi)
Conversion from Hz to Angular Hz.
Definition: constants.h:412
constexpr auto joule2kaycm(T x) -> decltype(x/kaycm2freq(h))
Conversion from Joule to cm-1.
Definition: constants.h:526
auto sind(T x) -> decltype(std::sin(deg2rad(x)))
Returns the sine of the deg2rad of the input.
Definition: constants.h:346
auto asind(T x) -> decltype(rad2deg(std::asin(x)))
Returns rad2deg of the arc-sine of the input.
Definition: constants.h:364
auto cosd(T x) -> decltype(std::cos(deg2rad(x)))
Returns the cosine of the deg2rad of the input.
Definition: constants.h:340
auto atan2d(T1 y, T2 x) -> decltype(rad2deg(std::atan2(y, x)))
Returns rad2deg of the arc-tangent of inputs #T1/#T2
Definition: constants.h:376
constexpr auto hz2joule(T x) -> decltype(x *h)
Conversion from MHz to Joule.
Definition: constants.h:502
constexpr auto kaycm2joule(T x) -> decltype(x *kaycm2freq(h))
Conversion from cm-1 to Joule.
Definition: constants.h:496
constexpr auto kelvin2celsius(T x) -> decltype(x - 273.15)
Conversion from K to C.
Definition: constants.h:466
auto acosd(T x) -> decltype(rad2deg(std::acos(x)))
Returns rad2deg of the arc-cosine of the input.
Definition: constants.h:358
constexpr auto freq2angcm(T x) -> decltype(x/kaycm2freq(inv_two_pi))
Conversion from Hz to Angular wavenumber.
Definition: constants.h:400
constexpr auto angcm2freq(T x) -> decltype(kaycm2freq(inv_two_pi))
Conversion from Angular wavenumber to Hz.
Definition: constants.h:394
constexpr auto freq2wavelen(T x) -> decltype(c/x)
Conversion from Hz to wavelength.
Definition: constants.h:424
constexpr auto mhz_per_torr2hz_per_pa(T x) -> decltype(x *pa2torr(1e6))
Conversion from MHz/Torr to Hz/Pa.
Definition: constants.h:454
constexpr auto pa2atm(T x) -> decltype(x/101 '325.0)
Conversion from Pa to Atm.
Definition: constants.h:436
constexpr auto mhz2joule(T x) -> decltype(hz2joule(x) *1e6)
Conversion from MHz to Joule.
Definition: constants.h:508
constexpr auto kaycm2freq(T x) -> decltype(x *(100 *c))
Conversion from Kayser wavenumber to Hz.
Definition: constants.h:382
constexpr auto angstrom2meter(T x) -> decltype(x *1e-10)
Conversion from Å to meter.
Definition: constants.h:532
constexpr auto kaycm_per_atm2hz_per_pa(T x) -> decltype(x *kaycm2freq(pa2atm(1)))
Conversion from cm-1 per atmosphere to Hz per Pascal.
Definition: constants.h:484
constexpr auto hz_per_msquared2kaycm_per_cmsquared(T x) -> decltype(x *freq2kaycm(1e4))
Conversion from Hz per molecule per m^2 to cm-1 per molecule per cm^2.
Definition: constants.h:478
constexpr auto kelvin2joule(T x) -> decltype(x *k)
Conversion from Kelvin to Joule.
Definition: constants.h:514
constexpr auto meter2angstrom(T x) -> decltype(x *1e10)
Conversion from meter to Å
Definition: constants.h:538
constexpr auto rad2deg(T x) -> decltype(x/one_degree_in_radians)
Converts radians to degrees.
Definition: constants.h:334
constexpr auto joule2hz(T x) -> decltype(x/h)
Conversion from Hz to Joule.
Definition: constants.h:520
constexpr auto pa2torr(T x) -> decltype(x/atm2pa(1.0/760.0))
Conversion from Pa to Torr.
Definition: constants.h:448
constexpr auto celsius2kelvin(T x) -> decltype(x+273.15)
Conversion from C to K.
Definition: constants.h:460
constexpr auto atm2pa(T x) -> decltype(x *101 '325.0)
Conversion from Atm to Pa.
Definition: constants.h:430
constexpr auto angfreq2freq(T x) -> decltype(x *inv_two_pi)
Conversion from Angular Hz to Hz.
Definition: constants.h:406
constexpr auto torr2pa(T x) -> decltype(x *atm2pa(1.0/760.0))
Conversion from Torr to Pa.
Definition: constants.h:442
constexpr auto kaycm_per_cmsquared2hz_per_msquared(T x) -> decltype(x *kaycm2freq(1e-4))
Conversion from cm-1 per molecule per cm^2 to Hz per molecule per m^2.
Definition: constants.h:472
constexpr auto deg2rad(T x) -> decltype(x *one_degree_in_radians)
Converts degrees to radians.
Definition: constants.h:328
constexpr auto freq2kaycm(T x) -> decltype(x/(100 *c))
Conversion from Hz to Kayser wavenumber.
Definition: constants.h:388
auto tand(T x) -> decltype(std::tan(deg2rad(x)))
Returns the tangent of the deg2rad of the input.
Definition: constants.h:352
constexpr auto wavelen2freq(T x) -> decltype(c/x)
Conversion from wavelength to Hz.
Definition: constants.h:418
auto atand(T x) -> decltype(rad2deg(std::atan(x)))
Returns rad2deg of the arc-tangent of the input.
Definition: constants.h:370
constexpr auto hz_per_pa2kaycm_per_atm(T x) -> decltype(x *freq2kaycm(atm2pa(1)))
Conversion from Hz per Pascal to cm-1 per atmosphere.
Definition: constants.h:490
X3 LineCenter QuadraticIndependent
Definition: constants.h:577
ENUMCLASS(TimeStep, char, hour, hours, h, minute, minutes, min, second, seconds, s) ENUMCLASS(LineShapeCoeff
Keep time options available to switch over them.
X3 LineStrength
Definition: constants.h:563
X3 LineCenter None
Definition: constants.h:576
#define u
#define v
#define w
#define c