![]() |
Equation of State Sub-Library: Version 0.910
|
00001 /* 00002 ------------------------------------------------------------------- 00003 00004 Copyright (C) 2006-2012, Andrew W. Steiner 00005 00006 This file is part of O2scl. 00007 00008 O2scl is free software; you can redistribute it and/or modify 00009 it under the terms of the GNU General Public License as published by 00010 the Free Software Foundation; either version 3 of the License, or 00011 (at your option) any later version. 00012 00013 O2scl is distributed in the hope that it will be useful, 00014 but WITHOUT ANY WARRANTY; without even the implied warranty of 00015 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00016 GNU General Public License for more details. 00017 00018 You should have received a copy of the GNU General Public License 00019 along with O2scl. If not, see <http://www.gnu.org/licenses/>. 00020 00021 ------------------------------------------------------------------- 00022 */ 00023 00024 #ifndef RMF_NUCLEUS_H 00025 #define RMF_NUCLEUS_H 00026 00027 #include <iostream> 00028 #include <string> 00029 #include <vector> 00030 #include <o2scl/uvector_tlate.h> 00031 #include <o2scl/umatrix_tlate.h> 00032 #include <o2scl/smart_interp.h> 00033 #include <o2scl/constants.h> 00034 #include <o2scl/part.h> 00035 #include <o2scl/rmf_eos.h> 00036 #include <o2scl/table_units.h> 00037 #include <o2scl/gsl_rkck.h> 00038 #include <o2scl/ode_funct.h> 00039 #include <o2scl/shared_ptr.h> 00040 00041 #ifndef DOXYGENP 00042 namespace o2scl { 00043 #endif 00044 00045 /** \brief Spherical closed-shell nuclei with a relativistic 00046 mean-field model in the Hartree approximation 00047 00048 This code is very experimental. 00049 00050 This class is based on a code developed by C.J. Horowitz and 00051 B.D. Serot, and used in \ref Horowitz81 which was then adapted 00052 by P.J. Ellis and used in \ref Heide94 and \ref Prakash94. Ellis 00053 and A.W. Steiner adapted it for the parameterization in in \ref 00054 rmf_eos for \ref Steiner05b, and then converted to C++ by 00055 Steiner afterwards. 00056 00057 The standard usage is something like: 00058 \code 00059 rmf_nucleus rn; 00060 o2scl_hdf::rmf_load(rn.rmf,"NL4"); 00061 rn.run_nucleus(82,208,0,0); 00062 cout << rn.rnrp << endl; 00063 \endcode 00064 which computes the structure of \f$ ^{208}\mathrm{Pb} \f$ and 00065 outputs the neutron skin thickness using the model \c 'NL4'. 00066 00067 Potential exceptions are 00068 - Failed to converge 00069 - Failed to solve meson field equations 00070 - Energy not finite (usually a problem in the equation of 00071 state) 00072 - Energy not finite in final calculation 00073 - Function \ref iterate() called before \ref init_run() 00074 - Not a closed-shell nucleus 00075 00076 The initial level pattern is 00077 \verbatim 00078 1 S 1/2 00079 // 2 nucleons 00080 1 P 3/2 00081 1 P 1/2 00082 // 8 nucleus 00083 1 D 5/2 00084 1 D 3/2 00085 2 S 1/2 00086 // 20 nucleons 00087 1 F 7/2 00088 // 28 nucleons 00089 1 F 5/2 00090 2 P 3/2 00091 2 P 1/2 00092 // 40 nucleons 00093 1 G 9/2 00094 // 50 nucleus 00095 1 G 7/2 00096 2 D 5/2 00097 1 H 11/2 00098 2 D 3/2 00099 3 S 1/2 00100 // 82 nucleons 00101 1 H 9/2 00102 2 F 7/2 00103 1 I 13/2 00104 2 F 5/2 00105 3 P 3/2 00106 3 P 1/2 00107 // 126 nucleons 00108 2 G 9/2 00109 1 I 11/2 00110 1 J 15/2 00111 3 D 5/2 00112 4 S 1/2 00113 2 G 7/2 00114 3 D 3/2 00115 // 184 nucleons 00116 \endverbatim 00117 00118 Below, \f$ \alpha \f$ is a generic index for the isospin, the 00119 radial quantum number \f$ n \f$ and the angular quantum numbers 00120 \f$ \kappa \f$ and \f$ m \f$. The meson fields are \f$ 00121 \sigma(r), \omega(r) \f$ and \f$ \rho(r) \f$. The baryon density 00122 is \f$ n(r) \f$, the neutron and proton densities are \f$ n_n(r) 00123 \f$ and \f$ n_p(r) \f$, and the baryon scalar density is \f$ 00124 n_s(r) \f$. 00125 The nucleon field equations are 00126 \f{eqnarray*} 00127 F^{\prime}_{\alpha}(r)- \frac{\kappa}{r} F_{\alpha}(r) 00128 + \left[ \varepsilon_{\alpha} - g_{\omega} \omega(r) 00129 - t_{\alpha} g_{\rho} \rho(r) - (t_{\alpha}+\frac{1}{2}) e A(r) - M 00130 + g_{\sigma} \sigma(r) \right] G_{\alpha}(r) &=& 0 \\ 00131 G^{\prime}_{\alpha}(r)+ \frac{\kappa}{r} G_{\alpha}(r) 00132 - \left[ \varepsilon_{\alpha} - g_{\omega} \omega(r) 00133 - t_{\alpha} g_{\rho} \rho(r) - (t_{\alpha}+\frac{1}{2}) e A(r) + M 00134 - g_{\sigma} \sigma(r) \right] F_{\alpha}(r) &=& 0 00135 \f} 00136 where \f$ t_{\alpha} \f$ is 1/2 for protons and -1/2 for 00137 neutrons. The meson field equations are 00138 \f{eqnarray*} 00139 \sigma^{\prime \prime}(r) + \frac{2}{r} \sigma^{\prime}(r) 00140 - m_{\sigma}^2 \sigma &=& - g_{\sigma} n_s(r) + b M g_{\sigma}^3 00141 \sigma^2 + c g_{\sigma}^4 \sigma^3 - g_{\rho}^2 \rho^2 00142 \frac{\partial f}{\partial \sigma} \\ 00143 \omega^{\prime \prime}(r) + \frac{2}{r} \omega^{\prime}(r) 00144 - m_{\omega}^2 \omega &=& - g_{\omega} n(r) + 00145 \frac{\zeta}{6} g_{\omega}^4 \omega^3 + g_{\rho}^2 \rho^2 00146 \frac{\partial f}{\partial \omega} \\ 00147 \rho^{\prime \prime}(r) + \frac{2}{r} \rho^{\prime}(r) 00148 - m_{\rho}^2 \rho &=& - \frac{g_{\rho}}{2} 00149 \left[n_n(r)-n_p(r)\right] + 2 g_{\rho}^2 \rho f + \frac{\xi}{6} 00150 g_{\rho}^4 \rho^3 00151 \f} 00152 and the Coulomb field equation is 00153 \f[ 00154 A^{\prime \prime}(r) + \frac{2}{r} A^{\prime}(r) = - e n_p(r) 00155 \f] 00156 00157 The densities are 00158 \f{eqnarray*} 00159 n_s &=& \sum_{\alpha} \left\{ \int d^3 r \left[ g(r)^2-f(r)^2 00160 \right] \right\} \\ 00161 n &=& \sum_{\alpha} \left\{ \int d^3 r \left[ g(r)^2+f(r)^2 00162 \right] \right\} \\ 00163 n_i &=& \sum_{\alpha} \left\{ t_{\alpha} \int d^3 r 00164 \left[ g(r)^2-f(r)^2 \right] \right\} \\ 00165 n_c &=& \sum_{\alpha} \left\{ \left[t_{\alpha}+\frac{1}{2}\right] 00166 \int d^3 r \left[ g(r)^2-f(r)^2 \right] \right\} 00167 \f} 00168 00169 Using the Green function 00170 \f[ 00171 D(r,r^{\prime},m_i) = \frac{-1}{m_i r r^{\prime}} 00172 \sinh (m_i r_{<}) \exp (-m_i r_{>}) 00173 \f] 00174 one can write the meson field 00175 \f[ 00176 \sigma(r) = \int_0^{\infty} r^{\prime 2}~d r^{\prime} 00177 \left[ - g_{\sigma} \rho_{s}(r) \right] 00178 D\left(r,r^{\prime},m_{\sigma}\right) 00179 \f] 00180 00181 When \f$ r > r^{\prime} \f$, and setting 00182 \f$ x = r^{\prime}-r \f$, the Green function is 00183 \f[ 00184 D = \frac{-1}{m_i r (x+r)} \sinh \left[m_i (x+r)\right] 00185 \exp (-m_i r) 00186 \f] 00187 When \f$ x >> r \f$, 00188 00189 00190 The total energy is 00191 \f[ 00192 E = \sum_{\alpha} \varepsilon_{\alpha} (2 j_{\alpha}+1) 00193 - \frac{1}{2} \int d^{3} r 00194 \left[- g_{\sigma} \sigma(r) \rho_s(r) + g_{\omega} \omega(r) 00195 \rho(r) + \frac{1}{2} g_{\rho} \rho(r) + e A(r) n_p(r) \right] 00196 \f] 00197 00198 The charge density is the proton density folded with the 00199 charge density of the proton, i.e. 00200 \f[ 00201 \rho_{\mathrm{ch}}(r) = \int d^{3} r^{\prime} 00202 \rho_{\mathrm{prot}}(r-r^{\prime}) \rho_p(r) 00203 \f] 00204 where the proton charge density is assumed to be of the form 00205 \f[ 00206 \rho_{\mathrm{prot}}(r) = \frac{\mu^3}{8 \pi} \exp \left( 00207 - \mu |r|\right) 00208 \f] 00209 and the parameter \f$ \mu = (0.71)^{1/2}~\mathrm{GeV} \f$ (see 00210 Eq. 20b in \ref Horowitz81). The default value of \ref a_proton 00211 is the value of \f$ \mu \f$ converted into \f$ \mathrm{fm}^{-1} 00212 \f$. 00213 00214 \todo Better documentation 00215 \todo Convert energies() to use EOS and possibly 00216 replace sigma_rhs() and related functions by the associated 00217 field equation method of rmf_eos. 00218 00219 \todo Document hw=3.923+23.265/cbrt(atot); 00220 00221 \comment 00222 the hw=blah blah term for the CM correction is discussed 00223 a bit in Negele, PRC 1 (1970) 1260, esp. see Eq. 2.30 but 00224 the numerical coefficients are different here. 00225 \endcomment 00226 00227 \future Sort energy levels at the end by energy 00228 \future Improve the numerical methods 00229 \future Make the neutron and proton orbitals more configurable 00230 \future Generalize to \f$ m_n \neq m_p \f$ . 00231 \future Allow more freedom in the integrations 00232 \future Consider converting everything to inverse fermis. 00233 \future Convert to zero-indexed arrays 00234 \future Warn when the level ordering is wrong, and unoccupied levels 00235 are lower energy than occupied levels 00236 */ 00237 class rmf_nucleus { 00238 00239 #ifndef DOXYGEN_INTERNAL 00240 00241 protected: 00242 00243 /// A convenient struct for the solution of the Dirac equations 00244 typedef struct { 00245 /// Eigenvalue 00246 double eigen; 00247 /// Quantum number \f$ \kappa \f$ . 00248 double kappa; 00249 /// The meson fields 00250 umatrix *fields; 00251 /// Desc 00252 uvector *varr; 00253 } odparms; 00254 00255 /// The total number of shells stored internally 00256 static const int n_internal_levels=29; 00257 00258 #endif 00259 00260 public: 00261 00262 rmf_nucleus(); 00263 00264 ~rmf_nucleus(); 00265 00266 /** \name Basic operation 00267 */ 00268 //@{ 00269 /// A shell of nucleons for \ref rmf_nucleus 00270 typedef struct { 00271 /// Degeneracy \f$ 2 j+1 \f$ . 00272 int twojp1; 00273 /// \f$ \kappa \f$ 00274 int kappa; 00275 /// Energy eigenvalue 00276 double energy; 00277 /// Isospin ( \f$ +1/2 \f$ or \f$ -1/2 \f$ . 00278 double isospin; 00279 /// Angular momentum-spin state \f$ ^{2s+1} \ell_{j} \f$ 00280 std::string state; 00281 /// Matching radius (in fm) 00282 double match_point; 00283 /// Desc. 00284 double eigen; 00285 /// Desc. 00286 double eigenc; 00287 /// Number of nodes in the wave function 00288 int nodes; 00289 } shell; 00290 00291 /** \brief Computes the structure of a nucleus with the specified 00292 number of levels 00293 00294 Note that \ref rmf must be set before run_nucleus() is called. 00295 00296 This calls init_run(), and then iterate() until \c iconverged is 00297 1, and then post_converge(). 00298 */ 00299 void run_nucleus(int nucleus_Z, int nucleus_N, int unocc_Z, int unocc_N); 00300 00301 /// Set output level 00302 void set_verbose(int v) { verbose=v; }; 00303 //@} 00304 00305 /** \name Lower-level interface 00306 */ 00307 //@{ 00308 /** \brief Initialize a run 00309 00310 Note that \ref rmf must be set before run_nucleus() is called. 00311 */ 00312 void init_run(int nucleus_Z, int nucleus_N, int unocc_Z, int unocc_N); 00313 00314 /// Perform an iteration 00315 void iterate(int nucleus_Z, int nucleus_N, int unocc_Z, int unocc_N, 00316 int &iconverged); 00317 00318 /// After convergence, make CM corrections, etc. 00319 int post_converge(int nucleus_Z, int nucleus_N, int unocc_Z, int unocc_N); 00320 00321 //@} 00322 00323 /// \name Results 00324 //@{ 00325 /** \brief Get the radial profiles 00326 00327 The profiles are calculated each iteration by iterate(). 00328 */ 00329 o2_shared_ptr<table_units>::type get_profiles() { return profiles; }; 00330 00331 /** \brief The finaal charge densities 00332 */ 00333 o2_shared_ptr<table_units>::type get_chden() { return chden_table; }; 00334 00335 /// The number of levels 00336 int nlevels; 00337 00338 /** \brief The levels (protons first, then neutrons) 00339 00340 An array of size \ref nlevels 00341 */ 00342 std::vector<shell> levels; 00343 00344 /// The number of unoccupied levels (equal to \c unocc_Z + \c unocc_N) 00345 int nuolevels; 00346 00347 /** \brief Information on the last convergence error 00348 */ 00349 int last_conv; 00350 00351 /** \brief The unoccupied levels (protons first, then neutrons) 00352 00353 An array of size \ref nuolevels 00354 */ 00355 std::vector<shell> unocc_levels; 00356 00357 /** \brief Surface tension (in \f$ \mathrm{fm}^{-3} \f$ ) 00358 00359 Computed in post_converge() or automatically in run_nucleus() 00360 */ 00361 double stens; 00362 00363 /** \brief Skin thickness (in fm) 00364 00365 Computed every iteration in iterate() 00366 or automatically in run_nucleus() 00367 */ 00368 double rnrp; 00369 00370 /** \brief Neutron RMS radius (in fm) 00371 00372 Computed every iteration in iterate() or automatically in 00373 run_nucleus() 00374 */ 00375 double rnrms; 00376 00377 /** \brief Proton RMS radius (in fm) 00378 00379 Computed every iteration in iterate() or automatically in 00380 run_nucleus() 00381 */ 00382 double rprms; 00383 00384 /** \brief Total energy (in MeV) 00385 00386 Computed every iteration in iterate() or automatically in 00387 run_nucleus() 00388 */ 00389 double etot; 00390 00391 /** \brief Charge radius (in fm) 00392 00393 Computed in post_converge() or automatically in run_nucleus() 00394 */ 00395 double r_charge; 00396 00397 /** \brief Charge radius corrected by the center of mass (in fm) 00398 00399 Computed in post_converge() or automatically in run_nucleus() 00400 */ 00401 double r_charge_cm; 00402 //@} 00403 00404 /** \name Equation of state 00405 */ 00406 //@{ 00407 /** \brief The default equation of state (default NL3) 00408 00409 This is set in the constructor to be the default 00410 model, NL3, using the function \ref load_nl3(). 00411 */ 00412 rmf_eos def_rmf; 00413 00414 /** \brief Set the base EOS to be used 00415 00416 The equation of state must be set before run_nucleus() or 00417 init_fun() are called, including the value of rmf_eos::mnuc. 00418 */ 00419 int set_eos(rmf_eos &r) { 00420 rmf=&r; 00421 return 0; 00422 } 00423 00424 /** \brief \ref thermo object for the EOS 00425 00426 This is just used as temporary storage. 00427 */ 00428 thermo hb; 00429 00430 /** \brief The neutron 00431 00432 The mass of the neutron is ignored and set by init_run() 00433 to be rmf_eos::mnuc from \ref rmf. 00434 */ 00435 fermion n; 00436 00437 /** \brief The proton 00438 00439 The mass of the proton is ignored and set by init_run() 00440 to be rmf_eos::mnuc from \ref rmf. 00441 */ 00442 fermion p; 00443 //@} 00444 00445 /** \name Numeric configuration 00446 */ 00447 //@{ 00448 /** \brief If true, call the error handler if the routine does not 00449 converge or reach the desired tolerance (default true) 00450 00451 If this is false, the function proceeds normally and 00452 may provide convergence information in \ref last_conv. 00453 */ 00454 bool err_nonconv; 00455 00456 /// The array type for the ODE solver 00457 typedef double arr_t[2]; 00458 00459 /// Set the stepper for the Dirac differential equation 00460 void set_step(ode_step<ode_funct<arr_t>,arr_t> &step) { 00461 ostep=&step; 00462 return; 00463 } 00464 00465 /// Maximum number of total iterations (default 70) 00466 int itmax; 00467 00468 /** \brief Maximum number of iterations for solving the meson 00469 field equations (default 10000) 00470 */ 00471 int meson_itmax; 00472 00473 /** \brief Maximum number of iterations for solving the Dirac 00474 equations (default 100) 00475 */ 00476 int dirac_itmax; 00477 00478 /// Tolerance for Dirac equations (default \f$ 5 \times 10^{-3} \f$ ). 00479 double dirac_tol; 00480 00481 /** \brief Second tolerance for Dirac equations 00482 (default \f$ 5 \times 10^{-4} \f$ ). 00483 */ 00484 double dirac_tol2; 00485 00486 /// Tolerance for meson field equations (default \f$ 10^{-6} \f$ ). 00487 double meson_tol; 00488 //@} 00489 00490 /** \brief Initial guess structure 00491 00492 The initial guess for the meson field profiles is 00493 a set of fermi-dirac functions, i.e. 00494 \f[ 00495 \sigma(r)=\mathrm{sigma0}/ 00496 [1+\exp((r-\mathrm{fermi\_radius})/\mathrm{fermi\_width})] 00497 \f] 00498 */ 00499 typedef struct { 00500 /// Scalar field at r=0 00501 double sigma0; 00502 /// Vector field at r=0 00503 double omega0; 00504 /// Isuvector field at r=0 00505 double rho0; 00506 /// Coulomb field at r=0 00507 double A0; 00508 /// The radius for which the fields are half their central value 00509 double fermi_radius; 00510 /// The "width" of the Fermi-Dirac function 00511 double fermi_width; 00512 } initial_guess; 00513 00514 /** \brief Parameters for initial guess 00515 00516 Default is {310,240,-6,25.9,6.85,0.6} 00517 */ 00518 initial_guess ig; 00519 00520 /** \brief If true, use the generic ODE solver instead of ... 00521 */ 00522 bool generic_ode; 00523 00524 #ifndef DOXYGEN_INTERNAL 00525 00526 protected: 00527 00528 /// The factor for the charge density of the proton (default 4.2707297) 00529 double a_proton; 00530 00531 /// Load the default model NL3 into the given \ref rmf_eos object 00532 int load_nl3(rmf_eos &r); 00533 00534 /// The base EOS 00535 rmf_eos *rmf; 00536 00537 /** \brief The radial profiles 00538 */ 00539 o2_shared_ptr<table_units>::type profiles; 00540 00541 /** \brief The final charge densities 00542 */ 00543 o2_shared_ptr<table_units>::type chden_table; 00544 00545 /** \brief A pointer to the current vector of levels 00546 (either \ref levels or \ref unocc_levels) 00547 */ 00548 std::vector<shell> *levp; 00549 00550 /** \brief Control output 00551 */ 00552 int verbose; 00553 00554 /// The starting neutron levels 00555 shell neutron_shells[n_internal_levels]; 00556 00557 /// The starting proton levels 00558 shell proton_shells[n_internal_levels]; 00559 00560 /// The grid size 00561 static const int grid_size=300; 00562 00563 /// The grid step size 00564 double step_size; 00565 00566 /// The nucleon mass (automatically set in init_fun()) 00567 double mnuc; 00568 00569 /** \name The meson fields and field equations (protected) 00570 */ 00571 //@{ 00572 /// Values of the fields from the last iteration 00573 umatrix field0; 00574 00575 /// The values of the fields 00576 umatrix fields; 00577 00578 /// The Green's functions inside 00579 umatrix gin; 00580 00581 /// The Green's functions outside 00582 umatrix gout; 00583 00584 /// Scalar density RHS 00585 double sigma_rhs(double sig, double ome, double rho); 00586 00587 /// Vector density RHS 00588 double omega_rhs(double sig, double ome, double rho); 00589 00590 /// Iso-vector density RHS 00591 double rho_rhs(double sig, double ome, double rho); 00592 00593 /// Calculate meson Green's functions 00594 void mesint(); 00595 00596 /// Calculate meson fields 00597 void meson(double ic); 00598 00599 /// Solve for the meson profiles 00600 void meson_solve(); 00601 00602 /** \brief The grid index corresponding to the nuclear surface 00603 (computed by init_run()) 00604 */ 00605 int surf_index; 00606 //@} 00607 00608 /** \name Density information (protected) 00609 */ 00610 //@{ 00611 /// The densities 00612 umatrix xrho; 00613 00614 /// The proton scalar density 00615 uvector xrhosp; 00616 00617 /// The scalar field RHS 00618 uvector xrhos; 00619 00620 /// The vector field RHS 00621 uvector xrhov; 00622 00623 /// The iso-vector field RHS 00624 uvector xrhor; 00625 00626 /// Charge density 00627 uvector chden1; 00628 00629 /// Charge density 00630 uvector chdenc; 00631 00632 /// Baryon density 00633 uvector arho; 00634 //@} 00635 00636 /// Energy profile 00637 uvector energy; 00638 00639 /// Initialize the meson fields, the densities, etc. 00640 void init_meson_density(); 00641 00642 /// Calculate the energy profile 00643 void energies(double xpro, double xnu, double e); 00644 00645 /// Compute the center of mass correction 00646 void center_mass_corr(double atot); 00647 00648 /** \name Calculating the form factor, etc. (protected) 00649 */ 00650 //@{ 00651 /// Fold in proton form factor 00652 void pfold(double x, double &xrhof); 00653 00654 /// Function representing proton form factor 00655 double xpform(double x, double xp, double a); 00656 00657 /// Perform integrations for form factor 00658 void gauss(double xmin, double xmax, double x, double &xi); 00659 00660 /// Desc. 00661 double xrhop(double x1); 00662 00663 /// Interpolation object 00664 smart_interp_vec<uvector_base,uvector_const_subvector, 00665 uvector,uvector_alloc> *gi; 00666 //@} 00667 00668 /** \name Solving the Dirac equations (protected) 00669 */ 00670 //@{ 00671 /** \brief Solve the Dirac equations 00672 00673 Solves the Dirac equation in from 12 fm to the match point and 00674 then out from .04 fm and adjusts eigenvalue with 00675 \f[ 00676 \Delta \varepsilon = -g(r=\mathrm{match\_point}) 00677 \times (f^{+}-f^{-}) 00678 \f] 00679 */ 00680 void dirac(int ilevel); 00681 00682 /// Take a step in the Dirac equations 00683 void dirac_step(double &x, double h, double eigen, 00684 double kappa, uvector &varr); 00685 00686 /// The form of the Dirac equations for the ODE stepper 00687 int odefun(double x, size_t nv, const arr_t &y, 00688 arr_t &dydx, odparms &op); 00689 00690 /// Compute the fields for the Dirac equations 00691 void field(double x, double &s, double &v, uvector_base &varr); 00692 00693 /// The default stepper 00694 gsl_rkck<ode_funct<arr_t>,arr_t,arr_t,array_alloc<arr_t> > def_step; 00695 00696 /// The ODE stepper 00697 ode_step<ode_funct<arr_t>,arr_t> *ostep; 00698 //@} 00699 00700 /// True if init() has been called 00701 bool init_called; 00702 00703 /** \brief Integrate the Dirac equations using a simple 00704 inline 4th order Runge-Kutta 00705 */ 00706 double gunt(double x, double g1, double f1, double &funt, 00707 double eigen, double kappa, uvector &varr); 00708 00709 /// ODE functions 00710 double ode_y[2]; 00711 00712 /// ODE derivatives 00713 double ode_dydx[2]; 00714 00715 /// ODE errors 00716 double ode_yerr[2]; 00717 00718 /// \name Gauss-Legendre integration points and weights 00719 //@{ 00720 double x12[6], w12[6]; 00721 double x100[50], w100[50]; 00722 //@} 00723 00724 #endif 00725 00726 }; 00727 00728 #ifndef DOXYGENP 00729 } 00730 #endif 00731 00732 #endif
Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).