Changeset 2034 for branches/c08_branch/source/hydro_bauman.cpp
- Timestamp:
- 05/10/08 09:03:02 (6 months ago)
- Files:
-
- 1 modified
-
branches/c08_branch/source/hydro_bauman.cpp (modified) (34 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/c08_branch/source/hydro_bauman.cpp
r1771 r2034 9 9 /* In this version, quantities that would normally cause a 64-bit floating point processor */ 10 10 /* to either underflow or overflow are evaluated using logs instead of floating point math. */ 11 /* This allows us to use an upper princip lequantum number `n' greater than which the */11 /* This allows us to use an upper principal quantum number `n' greater than which the */ 12 12 /* other version begins to fail. The trade-off is, of course, lower accuracy */ 13 13 /* ( or is it precision ). We use LOG_10 for convenience. */ … … 34 34 typedef struct t_mxq mxq; 35 35 36 /*lint -e790 integral to float */37 36 /************************************************************************************************/ 38 37 /* these routines were written by Robert Bauman */ … … 43 42 /* The recombination coefficient is obtained from the */ 44 43 /* photoionization cross-section (see Burgess 1965). */ 45 /* We have ,*/44 /* We have: */ 46 45 /* */ 47 46 /* - - l'=l+1 */ … … 124 123 /* r^2 */ 125 124 /* */ 126 /* Similiarly, the recombination co -efficient satisfies */125 /* Similiarly, the recombination coefficient satisfies */ 127 126 /************************************************************************************************/ 128 127 … … 144 143 /** \verbatim 145 144 ************************* for LOG version of the file *************************************** 146 In this version, quantit es that would normal cause a 64-bit floating point processor145 In this version, quantities that would normal cause a 64-bit floating point processor 147 146 to underflowed or overflow on intermediate values (ones internal to the calculation) 148 are evaluated using logs. This allows us to use an upper princip lequantum number `n'147 are evaluated using logs. This allows us to use an upper principal quantum number `n' 149 148 greater than 50 which is where the other version begins to fail. The trade-off is, 150 149 of course, lower accuracy( or is it precision ) and perhaps speed. … … 352 351 353 352 /** \verbatim 354 This routine, hri(), calculates the hydrogen radial inte rgral,353 This routine, hri(), calculates the hydrogen radial integral, 355 354 for the transition n,l --> n',l' 356 355 It is, of course, dimensionless. … … 373 372 374 373 /**\verbatim 375 This routine, hri_log10(), calculates the hydrogen radial inte rgral,374 This routine, hri_log10(), calculates the hydrogen radial integral, 376 375 for the transition n,l --> n',l' 377 376 It is, of course, dimensionless. … … 539 538 540 539 /************************************************************************/ 541 /* (4.0/3.0) * PI * FINE_STRUCTURE_CONSTANT * BO RHRADIUS * BORHRADIUS */540 /* (4.0/3.0) * PI * FINE_STRUCTURE_CONSTANT * BOHRRADIUS * BOHRRADIUS */ 542 541 /* */ 543 542 /* */ … … 742 741 /* K = k / Z ::: k^2 = ejected-electron-energy (Ryd) */ 743 742 double K, 744 /* princip lequantum number */743 /* principal quantum number */ 745 744 long int n, 746 745 /* angular momentum quantum number */ 747 746 long int l, 748 /* Temp erary storage for intermediate */747 /* Temporary storage for intermediate */ 749 748 /* results of the recursive routine */ 750 749 double *rcsvV … … 754 753 755 754 long int lp = 0; /* l' */ 756 double sigma = 0.; /* Sum in brocklehurst eq. 3.13 */755 double sigma = 0.; /* Sum in Brocklehurst eq. 3.13 */ 757 756 758 757 ASSERT( n > 0 ); … … 779 778 /* K = k / Z ::: k^2 = ejected-electron-energy (Ryd) */ 780 779 double K, 781 /* princip lequantum number */780 /* principal quantum number */ 782 781 long int n, 783 782 /* angular momentum quantum number */ 784 783 long int l, 785 /* Temp erary storage for intermediate */784 /* Temporary storage for intermediate */ 786 785 mxq *rcsvV_mxq 787 786 /* results of the recursive routine */ … … 791 790 792 791 long int lp = 0; /* l' */ 793 double sigma = 0.; /* Sum in brocklehurst eq. 3.13 */792 double sigma = 0.; /* Sum in Brocklehurst eq. 3.13 */ 794 793 795 794 ASSERT( n > 0 ); … … 835 834 long int l, 836 835 long int lp, 837 /* Temp erary storage for intermediate */836 /* Temporary storage for intermediate */ 838 837 /* results of the recursive routine */ 839 838 double *rcsvV … … 913 912 long int l, 914 913 long int lp, 915 /* Temp erary storage for intermediate */914 /* Temporary storage for intermediate */ 916 915 /* results of the recursive routine */ 917 916 mxq *rcsvV_mxq … … 1010 1009 long int l, 1011 1010 long int lp, 1012 /* Temp erary storage for intermediate */1011 /* Temporary storage for intermediate */ 1013 1012 /* results of the recursive routine */ 1014 1013 double *rcsvV … … 1104 1103 long int l, 1105 1104 long int lp, 1106 /* Temp erary storage for intermediate */1105 /* Temporary storage for intermediate */ 1107 1106 /* results of the recursive routine */ 1108 1107 mxq *rcsvV_mxq … … 1304 1303 /* */ 1305 1304 /* from the reference */ 1306 /* M. Broc Klehurst */1305 /* M. Brocklehurst */ 1307 1306 /* Mon. Note. R. astr. Soc. (1971) 153, 471-490 */ 1308 1307 /* */ … … 1320 1319 long int l, 1321 1320 long int lp, 1322 /* Temp erary storage for intermediate */1321 /* Temporary storage for intermediate */ 1323 1322 /* results of the recursive routine */ 1324 1323 double *rcsvV, … … 1471 1470 long int l, 1472 1471 long int lp, 1473 /* Temp erary storage for intermediate */1472 /* Temporary storage for intermediate */ 1474 1473 /* results of the recursive routine */ 1475 1474 mxq *rcsvV_mxq, … … 1680 1679 /* */ 1681 1680 /* from the reference */ 1682 /* M. Broc Klehurst */1681 /* M. Brocklehurst */ 1683 1682 /* Mon. Note. R. astr. Soc. (1971) 153, 471-490 */ 1684 1683 /* */ … … 1994 1993 /* */ 1995 1994 /* From reference; */ 1996 /* M. Broc Klehurst */1995 /* M. Brocklehurst */ 1997 1996 /* Mon. Note. R. astr. Soc. (1971) 153, 471-490 */ 1998 1997 /* */ … … 2022 2021 /* */ 2023 2022 /* */ 2024 /* the starting point for the recursion relations are ;*/2023 /* the starting point for the recursion relations are: */ 2025 2024 /* */ 2026 2025 /* */ 2027 /* equation : (3.18)*/2026 /* equation (3.18): */ 2028 2027 /* */ 2029 2028 /* | pi |(1/2) 8n */ … … 2031 2030 /* | 2 | (2n-1)! */ 2032 2031 /* */ 2033 /* equation : (3.20)*/2032 /* equation (3.20): */ 2034 2033 /* */ 2035 2034 /* exp(2n-2/K tan^(-1)(n K) */ … … 2039 2038 /* */ 2040 2039 /* */ 2041 /* equation : (3.20)*/2040 /* equation (3.20): */ 2042 2041 /* G(n,n-2; K,n-1) = (2n-2)(1+(n K)^2) n G(n,n-1; K,n) */ 2043 2042 /* */ 2044 2043 /* */ 2045 /* equation : (3.21)*/2044 /* equation (3.21): */ 2046 2045 /* */ 2047 2046 /* (1+(n K)^2) */ … … 2055 2054 long int l, 2056 2055 long int lp, 2057 /* Temp erary storage for intermediate */2056 /* Temporary storage for intermediate */ 2058 2057 /* results of the recursive routine */ 2059 2058 double *rcsvV … … 2115 2114 long int l, 2116 2115 long int lp, 2117 /* Temp erary storage for intermediate */2116 /* Temporary storage for intermediate */ 2118 2117 /* results of the recursive routine */ 2119 2118 mxq *rcsvV_mxq … … 2384 2383 if( n >= 70 ) 2385 2384 { 2386 fprintf(ioQQQ,"Princip leQuantum Number `n' too large.\n");2385 fprintf(ioQQQ,"Principal Quantum Number `n' too large.\n"); 2387 2386 cdEXIT(EXIT_FAILURE); 2388 2387 } … … 2399 2398 if( n <= np ) 2400 2399 { 2401 fprintf(ioQQQ," The princip lequantum numbers are such that n <= n'.\n");2400 fprintf(ioQQQ," The principal quantum numbers are such that n <= n'.\n"); 2402 2401 cdEXIT(EXIT_FAILURE); 2403 2402 } … … 2440 2439 if( n <= np ) 2441 2440 { 2442 fprintf(ioQQQ," The princip lequantum numbers are such that n <= n'.\n");2441 fprintf(ioQQQ," The principal quantum numbers are such that n <= n'.\n"); 2443 2442 cdEXIT(EXIT_FAILURE); 2444 2443 } … … 2521 2520 if( n <= nprime ) 2522 2521 { 2523 fprintf(ioQQQ," The princip lequantum numbers are such that n <= n'.\n");2522 fprintf(ioQQQ," The principal quantum numbers are such that n <= n'.\n"); 2524 2523 cdEXIT(EXIT_FAILURE); 2525 2524 } … … 2530 2529 /************************************************************************/ 2531 2530 /* hri() */ 2532 /* Calculate the hydrogen radial wavefunction inte rgral*/2531 /* Calculate the hydrogen radial wavefunction integral */ 2533 2532 /* for the dipole transition l'=l-1 or l'=l+1 */ 2534 2533 /* for the higher energy state n,l to the lower energy state n',l' */ … … 2551 2550 2552 2551 /************************************************************************/ 2553 /* This routine, hri(), calculates the hydrogen radial inte rgral,*/2552 /* This routine, hri(), calculates the hydrogen radial integral, */ 2554 2553 /* for the transition n,l --> n',l' */ 2555 2554 /* It is, of course, dimensionless. */ … … 2628 2627 /************************************************************************/ 2629 2628 /* hri_log10() */ 2630 /* Calculate the hydrogen radial wavefunction inte rgral*/2629 /* Calculate the hydrogen radial wavefunction integral */ 2631 2630 /* for the dipole transition l'=l-1 or l'=l+1 */ 2632 2631 /* for the higher energy state n,l to the lower energy state n',l' */ … … 4445 4444 return partsum; 4446 4445 } 4447 /*lint +e790 integral to float */
