TFluka/Fdblprc.h

   1 #ifndef FDBLPRC_H
   2 #define FDBLPRC_H 1
   3
   4 #include "Rtypes.h"
   5 #include "cfortran.h"
   6 extern "C" {
   7 //*$ create dblprc.add
   8 //*copy dblprc
   9 //*                                                                     *
  10 //*=== dblprc ==========================================================*
  11 //*                                                                     *
  12 //*---------------------------------------------------------------------*
  13 //*                                                                     *
  14 //*      dblprc: included in any routine, machine, mathematical and     *
  15 //*              physical constants plus global declarations            *
  16 //*                                                                     *
  17 //*  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  *
  18 //*  !!!! o n   m a c h i n e s   w h e r e   t h e   d o u b l e !!!!  *
  19 //*  !!!! p r e c i s i o n   i s   n o t   r e q u i r e d  r e -!!!!  *
  20 //*  !!!! m o v e   t h e   d o u b l e   p r e c i s i o n       !!!!  *
  21 //*  !!!! s t a t e m e n t,  s e t   k a l g n m = 1   a n d     !!!!  *
  22 //*  !!!! c h a n g e   a l l   n u m e r i c a l   c o n s -     !!!!  *
  23 //*  !!!! t a n t s   t o   s i n g l e   p r e c i s i o n       !!!!  *
  24 //*  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  *
  25 //*                                                                     *
  26 //*         kalgnm = real address alignment, 2 for double precision,    *
  27 //*                  1 for single precision                             *
  28 //*         kalch8 = character*8 address alignment wrt the precision    *
  29 //*                  defined by kalgnm (mostly 1 in all situations)     *
  30 //*         i2algn = integer*2 address alignment wrt the normal integer *
  31 //*                  precision (mostly 2, 4 for 64 bit integers)        *
  32 //*         anglgb = this parameter should be set equal to the machine  *
  33 //*                  "zero" with respect to unit                        *
  34 //*         anglsq = this parameter should be set equal to the square   *
  35 //*                  of anglgb                                          *
  36 //*         axcssv = this parameter should be set equal to the number   *
  37 //*                  for which unity is negligible for the machine      *
  38 //*                  accuracy                                           *
  39 //*         andrfl = "underflow" of the machine for floating point      *
  40 //*                  operation                                          *
  41 //*         avrflw = "overflow"  of the machine for floating point      *
  42 //*                  operation                                          *
  43 //*         ainfnt = code "infinite"                                    *
  44 //*         azrzrz = code "zero"                                        *
  45 //*         einfnt = natural logarithm of the code "infinite"           *
  46 //*         ezrzrz = natural logarithm of the code "zero"               *
  47 //*         excssv = natural logarithm of the code number for which     *
  48 //*                  unit is negligible                                 *
  49 //*         englgb = natural logarithm of the code "zero" with respect  *
  50 //*                  to unit                                            *
  51 //*         onemns = 1- of the machine, it is 1 - 2 x anglgb            *
  52 //*         onepls = 1+ of the machine, it is 1 + 2 x anglgb            *
  53 //*         csnnrm = maximum tolerable error on cosine normalization,   *
  54 //*                  u**2+v**2+w**2: assuming a typical anglgb relative *
  55 //*                  error on each component we would get 2xanglgb: use *
  56 //*                  4xanglgb to avoid too many normalizations          *
  57 //*         dmxtrn = "infinite" distance for transport (cm)             *
  58 //*         rhflmn = minimal density for fluka (g/cm^3)                 *
  59 //*                                                                     *
  60 //*   "global" declarations:                                            *
  61 //*         lfluka = set to true for a real (full) fluka run            *
  62 //*         lgbias = set to true for a fully biased run                 *
  63 //*         lgbana = set to true for a fully analogue run               *
  64 //*         lflgeo = set to true when using the standard fluka geometry *
  65 //*         loflts = set to true for special off-line testing of speci- *
  66 //*                  fic routines                                       *
  67 //*         lusrin = set to true if the user dependent initialization   *
  68 //*                  routine usrini has been called at least once       *
  69 //*         lnmgeo = set to true for a name-base geometry input         *
  70 //*         lnminp = set to true for a name-base fluka input            *
  71 //*         Lfrfmt = set to true for a free-format based Fluka input    *
  72 //*         lfdrtr = set to true for going in/out feeder/flukam at each *
  73 //*                  event                                              *
  74 //*                                                                     *
  75 //*---------------------------------------------------------------------*
  76 //*                                                                     *
  77 const Int_t kalgnm = 2;
  78 const Int_t kalch8 = 1;
  79 const Int_t i2algn = 2;
  80 const Double_t anglgb = 5.0e-16;
  81 const Double_t anglsq = 2.5e-31;
  82 const Double_t axcssv = 0.2e+16;
  83 const Double_t andrfl = 1.0e-38;
  84 const Double_t avrflw = 1.0e+38;
  85 const Double_t ainfnt = 1.0e+30;
  86 const Double_t azrzrz = 1.0e-30;
  87 const Double_t einfnt = +69.07755278982137e+00;
  88 const Double_t ezrzrz = -69.07755278982137e+00;
  89 const Double_t excssv = +35.23192357547063e+00;
  90 const Double_t englgb = -35.23192357547063e+00;
  91 const Double_t onemns = 0.999999999999999e+00;
  92 const Double_t onepls = 1.000000000000001e+00;
  93 const Double_t csnnrm = 2.0e-15;
  94 const Double_t dmxtrn = 1.0e+08;
  95 const Double_t rhflmn = 1.0e-10;
  96 //*
  97 //*======================================================================*
  98 //*======================================================================*
  99 //*=========                                                   ==========*
 100 //*=========    m a t h e m a t i c a l   c o n s t a n t s    ==========*
 101 //*=========                                                   ==========*
 102 //*======================================================================*
 103 //*======================================================================*
 104 //*                                                                      *
 105 //*   numerical constants (single precision):                            *
 106 //*                                                                      *
 107 //*         zersng = 0                                                   *
 108 //*                                                                      *
 109 //*   numerical constants (double precision):                            *
 110 //*                                                                      *
 111 //*         zerzer = 0                                                   *
 112 //*         oneone = 1                                                   *
 113 //*         twotwo = 2                                                   *
 114 //*         thrthr = 3                                                   *
 115 //*         foufou = 4                                                   *
 116 //*         fivfiv = 5                                                   *
 117 //*         sixsix = 6                                                   *
 118 //*         sevsev = 7                                                   *
 119 //*         eigeig = 8                                                   *
 120 //*         aninen = 9                                                   *
 121 //*         tenten = 10                                                  *
 122 //*         eleven = 11                                                  *
 123 //*         twelve = 12                                                  *
 124 //*         fiften = 15                                                  *
 125 //*         sixten = 16                                                  *
 126 //*         hlfhlf = 1/2                                                 *
 127 //*         onethi = 1/3                                                 *
 128 //*         onefou = 1/4                                                 *
 129 //*         onefiv = 1/5                                                 *
 130 //*         onesix = 1/6                                                 *
 131 //*         onesev = 1/7                                                 *
 132 //*         oneeig = 1/8                                                 *
 133 //*         twothi = 2/3                                                 *
 134 //*         thrfou = 3/4                                                 *
 135 //*         thrtwo = 3/2                                                 *
 136 //*         two2o3 = 2^2/3                                               *
 137 //*         pipipi = circumference / diameter                            *
 138 //*         twopip = 2 x pipipi                                          *
 139 //*         pip5o2 = 5/2 x pipipi                                        *
 140 //*         pipisq = pipipi x pipipi                                     *
 141 //*         pihalf = 1/2 x pipipi                                        *
 142 //*         erfa00 = erf (oo) = 1/2 x square root of pi                  *
 143 //*         sqtwpi = square root of 2xpi                                 *
 144 //*         eulero = eulero's constant                                   *
 145 //*         eulexp = exp ( eulero )                                      *
 146 //*         e1m2eu = exp ( 1 - 2 eulero )                                *
 147 //*         eneper = "e", base of natural logarithm                      *
 148 //*         sqrent = square root of "e"                                  *
 149 //*         sqrtwo = square root of  2                                   *
 150 //*         sqrthr = square root of  3                                   *
 151 //*         sqrfiv = square root of  5                                   *
 152 //*         sqrsix = square root of  6                                   *
 153 //*         sqrsev = square root of  7                                   *
 154 //*         sqrt12 = square root of 12                                   *
 155 //*         s2fwhm = 2 x square root of 2 x logarithm of 2               *
 156 //*                                                                      *
 157 //*----------------------------------------------------------------------*
 158 //*
 159     const Float_t  zersng = 0.e+00;
 160     const Double_t zerzer = 0.e+00;
 161     const Double_t oneone = 1.e+00;
 162     const Double_t twotwo = 2.e+00;
 163     const Double_t thrthr = 3.e+00;
 164     const Double_t foufou = 4.e+00;
 165     const Double_t fivfiv = 5.e+00;
 166     const Double_t sixsix = 6.e+00;
 167     const Double_t sevsev = 7.e+00;
 168     const Double_t eigeig = 8.e+00;
 169     const Double_t aninen = 9.e+00;
 170     const Double_t tenten = 10.e+00;
 171     const Double_t eleven = 11.e+00;
 172     const Double_t twelve = 12.e+00;
 173     const Double_t fiften = 15.e+00;
 174     const Double_t sixten = 16.e+00;
 175     const Double_t hlfhlf = 0.5e+00;
 176     const Double_t onethi = oneone/thrthr;
 177     const Double_t onefou = oneone/foufou;
 178     const Double_t onefiv = oneone/fivfiv;
 179     const Double_t onesix = oneone/sixsix;
 180     const Double_t onesev = oneone/sevsev;
 181     const Double_t oneeig = oneone/eigeig;
 182     const Double_t twothi = twotwo/thrthr;
 183     const Double_t thrfou = thrthr/foufou;
 184     const Double_t thrtwo = thrthr/twotwo;
 185     const Double_t fouthr = foufou/thrthr;
 186     const Double_t pipipi = 3.141592653589793238462643383279e+00;
 187     const Double_t two2o3 = 1.587401051968199e+00; // 2006.3
 188     const Double_t twopip = 6.283185307179586476925286766559e+00;
 189     const Double_t pip5o2 = 7.853981633974483096156608458199e+00;
 190     const Double_t pipisq = 9.869604401089358618834490999876e+00;
 191     const Double_t pihalf = 1.570796326794896619231321691640e+00;
 192     const Double_t erfa00 = 0.886226925452758013649083741671e+00;
 193     const Double_t sqrtpi = 1.772453850905516027298167483341e+00;
 194     const Double_t sqtwpi = 2.506628274631000502415765284811e+00;
 195     const Double_t eulero = 0.577215664901532860606512e+00;
 196     const Double_t eulexp = 1.781072417990197985236504e+00;
 197     const Double_t eullog = -0.5495393129816448223376619e+00;
 198     const Double_t e1m2eu = 0.8569023337737540831433017e+00;
 199     const Double_t eneper = 2.718281828459045235360287471353e+00;
 200     const Double_t sqrent = 1.648721270700128146848650787814e+00;
 201     const Double_t sqrtwo = 1.414213562373095048801688724210e+00;
 202     const Double_t sqrthr = 1.732050807568877293527446341506e+00;
 203     const Double_t sqrfiv = 2.236067977499789696409173668731e+00;
 204     const Double_t sqrsix = 2.449489742783178098197284074706e+00;
 205     const Double_t sqrsev = 2.645751311064590590501615753639e+00;
 206     const Double_t sqrt12 = 3.464101615137754587054892683012e+00;
 207     const Double_t s2fwhm = 2.354820045030949e+00;
 208     const Double_t twolog = 0.693147180559945309417232121458e+00;
 209 //*
 210 //*======================================================================*
 211 //*======================================================================*
 212 //*=========                                                   ==========*
 213 //*=========       p h y s i c a l   c o n s t a n t s         ==========*
 214 //*=========                                                   ==========*
 215 //*======================================================================*
 216 //*======================================================================*
 217 //*                                                                      *
 218 //*   primary constants:                                                 *
 219 //*                                                                      *
 220 //*         clight = speed of light in cm s-1                            *
 221 //*         avogad = avogadro number                                     *
 222 //*         boltzm = k boltzmann constant (j k-1)                        *
 223 //*         amelgr = electron mass (g)                                   *
 224 //*         plckbr = reduced planck constant (erg s)                     *
 225 //*         elccgs = elementary charge (cgs unit)                        *
 226 //*         elcmks = elementary charge (mks unit)                        *
 227 //*         epsil0 = vacuum dielectric constant (MKS unit)               *
 228 //*         amugrm = atomic mass unit (g)                                *
 229 //*         ammumu = muon    mass (amu)                                  *
 230 //*         amprmu = proton  mass (amu)                                  *
 231 //*         amnemu = neutron mass (amu)                                  *
 232 //*                                                                      * //*          Note:                                                       *
 233 //*          e[MKS]^2/ (4 pi epsilon0 hbar[MKS] c[MKS]) = alpha = 1/137  *
 234 //*          e[CGS]^2/ (hbar[CGS] c[CGS]) = alpha = 1/137                *
 235 //*          with c[MKs] = 10^-2 c[CGS], hbar[MKS] = 10^-7 hbar[CGS]     *
 236 //*   derived constants:                                                 *
 237 //*                                                                      *
 238 //*         alpfsc = fine structure constant  = e^2/(hbar c) (cgs units) *
 239 //*         amelct = electron mass (gev) = 10^-16amelgr clight^2 / elcmks*
 240 //*         amugev = atomic mass unit (gev) = 10^-16amugrm clight^2      *
 241 //*                                           / elcmks                   *
 242 //*         ammuon = muon    mass (gev) = ammumu * amugev                *
 243 //*         amprtn = proton  mass (gev) = amprmu * amugev                *
 244 //*         amntrn = neutron mass (gev) = amnemu * amugev                *
 245 //*         amdeut = deuteron mass (gev)                                 *
 246 //*         amalph = alpha    mass (gev) (derived from the excess mass   *
 247 //*                  and an (approximate) atomic binding not a really    *
 248 //*                  measured constant)                                  *
 249 //*         cougfm = e^2 (gev fm) = elccgs^2 / elcmks * 10^-7 * 10^-9    *
 250 //*                * 10^13 (10^..=erg cm->joule cm->gev cm->gev fm       *
 251 //*                it is equal to 0.00144 gev fm                         *
 252 //*         fscto2 = (fine structure constant)^2                         *
 253 //*         fscto3 = (fine structure constant)^3                         *
 254 //*         fscto4 = (fine structure constant)^4                         *
 255 //*         plabrc = reduced planck constant times the light velocity    *
 256 //*                  expressed in gev fm                                 *
 257 //*         rclsel = classical electron radius (cm) = e^2 / (m_e c^2)    *
 258 //*         bltzmn = k boltzmann constant in gev k-1                     *
 259 //*         a0bohr = bohr radius, hbar^2 / ( m_e e^2) (fm) = plabrc**2   *
 260 //*                / amelct / cougfm, or equivalently,                   *
 261 //*                plabrc / alpfsc / amelct                              *
 262 //*         gfohb3 = fermi constant, g_f/(hbar c)^3, in gev^-2           *
 263 //*         gfermi = fermi constant in gev fm^3                          *
 264 //*         sin2tw = sin^2 theta_weinberg                                *
 265 //*         prmgnm = proton  magnetic moment (magneton)                  *
 266 //*         anmgnm = neutron magnetic moment (magneton)                  *
 267 //*         s0thms = sigma_0 Thomson, 8/3 pi r_e^2 (mb)                  *
 268 //*                                                                      *
 269 //*   astronomical constants:                                            *
 270 //*                                                                      *
 271 //*         rearth = earth equatorial radius (cm)                        *
 272 //*         auastu = astronomical unit       (cm)                        *
 273 //*                                                                      *
 274 //*   conversion constants:                                              *
 275 //*                                                                      *
 276 //*         gevmev = from gev to mev                                     *
 277 //*         emvgev = from mev to gev                                     *
 278 //*         gev2ev = from gev to  ev                                     *
 279 //*         ev2gev = from ev  to gev                                     *
 280 //*         algvmv = from gev to mev, log                                *
 281 //*         raddeg = from radians to degrees                             *
 282 //*         degrad = from degrees to radians                             *
 283 //*         gevomg = from (photon) energy [gev] in 2pi x frequency [s^-1]*
 284 //*         cmq2mb = from square centimetres to millibarns               *
 285 //*                                                                      *
 286 //*   useful constants:                                                  *
 287 //*                                                                      *
 288 //*         fertho = constant to be used in the fermi-thomas approxima-  *
 289 //*                  ted expression for atomic binding energies          *
 290 //*         expebn = exponent to be used in the fermi-thomas approxima-  *
 291 //*                  ted expression for atomic binding energies          *
 292 //*                    b_atomic (z) = fertho x z^expebn (gev)            *
 293 //*         bexc12 = fermi-thomas approximated expression for 12-c ato-  *
 294 //*                  mic binding energies (gev)                          *
 295 //*         amunmu = difference between the atomic and nuclear mass units*
 296 //*         amuc12 = "nuclear" mass unit = 1/12 m_nucl (12-c),           *
 297 //*                  m_nucl (12-c) = m_atom (12-c) - 6 m_e + b_atom(12-c)*
 298 //*                                                                      *
 299 //*----------------------------------------------------------------------*
 300 //*
 301     const Double_t clight = 2.99792458e+10;
 302     const Double_t avogad = 6.0221367e+23;
 303     const Double_t boltzm = 1.380658e-23;
 304     const Double_t amelgr = 9.1093897e-28;
 305     const Double_t plckbr = 1.05457266e-27;
 306     const Double_t elccgs = 4.8032068e-10;
 307     const Double_t elcmks = 1.60217733e-19;
 308     const Double_t epsil0 = 8.854187817e-12; // 2006.3
 309     const Double_t amugrm = 1.6605402e-24;
 310     const Double_t ammumu = 0.113428913e+00;
 311     const Double_t amprmu = 1.007276470e+00;
 312     const Double_t amnemu = 1.008664904e+00;
 313 //* const Double_t alpfsc = 1.e+00 / 137.035989561e+00
 314 //* const Double_t fscto2 = alpfsc * alpfsc
 315 //* const Double_t fscto3 = fscto2 * alpfsc
 316 //* const Double_t fscto4 = fscto3 * alpfsc
 317 //*    it is important to set the electron mass exactly with the same
 318 //*    rounding as in the mass tables, so use the explicit expression
 319 //* const Double_t amelct = 1.e-16 * amelgr * clight * clight / elcmks
 320 //*    it is important to set the amu mass exactly with the same
 321 //*    rounding as in the mass tables, so use the explicit expression
 322 //* const Double_t amugev = 1.e-16 * amugrm * clight * clight / elcmks
 323 //*    it is important to set the muon,proton,neutron masses exactly with
 324 //*    the same rounding as in the mass tables, so use the explicit
 325 //*    expression
 326 //* const Double_t ammuon = ammumu * amugev
 327 //* const Double_t amprtn = amprmu * amugev
 328 //* const Double_t amntrn = amnemu * amugev
 329 //* const Double_t rclsel = elccgs * elccgs / clight / clight / amelgr
 330 //* const Double_t bltzmn = boltzm / elcmks * 1.e-09
 331 const Double_t alpfsc = 7.2973530791728595e-3;
 332 const Double_t fscto2 = 5.3251361962113614e-5;
 333 const Double_t fscto3 = 3.8859399018437826e-7;
 334 const Double_t fscto4 = 2.8357075508200407e-9;
 335 const Double_t plabrc = 0.197327053e+00;
 336 const Double_t amelct = 0.51099906e-3;
 337 const Double_t amugev = 0.93149432e+00;
 338 const Double_t ammuon = 0.105658389e+00;
 339 const Double_t amprtn = 0.93827231e+00;
 340 const Double_t amntrn = 0.93956563e+00;
 341 const Double_t amdeut = 1.87561339e+00;
 342 const Double_t amalph = 3.72738025692891e+00;
 343 const Double_t cougfm = elccgs*elccgs/elcmks*(1.e-7)*(1.e+13)*(1.e-9);
 344 const Double_t rclsel = 2.8179409183694872e-13;
 345 const Double_t alamb0 = twotwo * pipipi * rclsel / alpfsc;
 346 const Double_t bltzmn = 8.617385e-14;
 347 const Double_t a0bohr = plabrc/alpfsc/amelct;
 348 const Double_t gfohb3 = 1.16639e-5;
 349 const Double_t gfermi = gfohb3*plabrc*plabrc*plabrc;
 350 const Double_t sin2tw = 0.2319e+00;
 351 const Double_t prmgnm = 2.792847386e+00;
 352 const Double_t anmgnm = -1.91304275e+00;
 353 const Double_t rearth = 6.378140e+8;
 354 const Double_t auastu = 1.4959787066e+13;
 355 const Double_t gevmev = 1.0e+3;
 356 const Double_t ev2gev = 1.0e-9;
 357 const Double_t gev2ev = 1.0e+9;
 358 const Double_t emvgev = 1.0e-3;
 359 const Double_t cmq2mb = 1.0e+27;
 360 const Double_t fmb2ba = 1.0e-3;
 361 const Double_t bar2mb = 1.0e+3;
 362 const Double_t fmb2fs = 1.0e-1;
 363 const Double_t fms2mb = 1.0e+1;
 364 const Double_t algvmv = 6.90775527898214e+00;
 365 const Double_t raddeg = (180.e+00)/pipipi;
 366 const Double_t degrad = pipipi/(180.e+00);
 367 const Double_t gevomg = clight*(1.e+13)/plabrc;
 368 const Double_t s0thms = eigeig / thrthr * pipipi * rclsel * rclsel * cmq2mb;
 369 //*  old Fermi-Thomas parametrization of atomic binding energies:
 370 //*     const Double_t fertho = 15.73       e-9
 371 //*     const Double_t expebn = 7.e+00 / 3.e+00
 372 //*     const Double_t bexc12 = fertho * 65.41634134195703e+00
 373 //*  new Fermi-Thomas parametrization of atomic binding energies:
 374 const Double_t fertho = 14.33e-9;
 375 const Double_t expebn = 2.39e+00;
 376 const Double_t bexc12 = fertho*72.40715579499394e+00;
 377 const Double_t amunmu = hlfhlf*amelct-bexc12/12.e+00;
 378 const Double_t amuc12 = amugev-amunmu;
 379 //*  Old MeV units:
 380 const Double_t amemev = gevmev * amelct;
 381 //*
 382
 383 typedef struct {
 384     Int_t    lfluka;
 385     Int_t    lgbias;
 386     Int_t    lgbana;
 387     Int_t    lflgeo;
 388     Int_t    loflts;
 389     Int_t    lusrin;
 390     Int_t    lusrgl; // 2006.3
 391     Int_t    lnmgeo;
 392     Int_t    lnminp;
 393     Int_t    lfrfmt;
 394     Int_t    lfdrtr;
 395     Int_t    kflgeo;
 396     Int_t    kfldnr;
 397 } globalCommon;
 398 #define GLOBAL COMMON_BLOCK(GLOBAL,global)
 399 COMMON_BLOCK_DEF(globalCommon,GLOBAL);
 400 }
 401
 402 #endif