00001 /* 00002 ------------------------------------------------------------------- 00003 00004 Copyright (C) 2006, 2007, 2008, 2009, 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 #ifndef O2SCL_FERMION_H 00024 #define O2SCL_FERMION_H 00025 00026 #include <string> 00027 #include <iostream> 00028 #include <fstream> 00029 #include <cmath> 00030 00031 // For gsl_sf_fermi_dirac_int() 00032 #include <gsl/gsl_specfunc.h> 00033 00034 #include <o2scl/constants.h> 00035 #include <o2scl/funct.h> 00036 #include <o2scl/root.h> 00037 #include <o2scl/cern_mroot_root.h> 00038 #include <o2scl/part.h> 00039 00040 #ifndef DOXYGENP 00041 namespace o2scl { 00042 #endif 00043 00044 /** 00045 \brief Fermion class 00046 00047 This is a base class for the computation of fermionic statistics 00048 at zero temperature. The more general case of finite temperature 00049 is taken care of by \ref fermion_T objects. The primary 00050 functions are calc_mu_zerot() and calc_density_zerot() which 00051 compute all the thermodynamic quantities as a function of the 00052 chemical potential, or the density, respectively. 00053 00054 This class also adds two member data variables, \ref kf and \ref 00055 del, for the Fermi momentum and the gap. 00056 00057 \future Use hypot() and other more accurate functions for the 00058 analytic expressions for the zero temperature integrals. [Progress 00059 has been made, but there are probably other functions which may 00060 break down for small but finite masses and temperatures] 00061 00062 */ 00063 class fermion : virtual public part { 00064 00065 public: 00066 00067 /// Fermi momentum 00068 double kf; 00069 /// Gap 00070 double del; 00071 00072 /// Create a fermion with mass \c mass and degeneracy \c dof. 00073 fermion(double mass=0, double dof=0); 00074 00075 virtual ~fermion() { 00076 } 00077 00078 friend class io_tlate<fermion>; 00079 00080 /// \name Zero-temperature fermions 00081 //@{ 00082 /** 00083 \brief Calculate the Fermi momentum from the density 00084 00085 Uses the relation \f$ k_F = ( 6 \pi^2 n /g )^{1/3} \f$ 00086 */ 00087 int kf_from_density(); 00088 00089 /** 00090 \brief Energy density at T=0 from \ref kf and \ref ms 00091 00092 Calculates the integral 00093 \f[ 00094 \varepsilon = \frac{g}{2 \pi^2} \int_0^{k_F} k^2 00095 \sqrt{k^2+m^{* 2}} d k 00096 \f] 00097 */ 00098 int energy_density_zerot(); 00099 00100 /** 00101 \brief Pressure at T=0 from \ref kf and \ref ms 00102 00103 Calculates the integral 00104 \f[ 00105 P=\frac{g}{6 \pi^2} \int_0^{k_F} \frac{k^4}{\sqrt{k^2+m^{* 2}}} d k 00106 \f] 00107 */ 00108 int pressure_zerot(); 00109 00110 /** 00111 \brief Zero temperature fermions from nu and ms 00112 00113 This function always returns \c gsl_success. 00114 */ 00115 virtual int calc_mu_zerot(); 00116 00117 /** 00118 \brief Zero temperature fermions from n and ms 00119 00120 This function always returns \c gsl_success. 00121 */ 00122 virtual int calc_density_zerot(); 00123 //@} 00124 00125 /// Return string denoting type ("fermion") 00126 virtual const char *type() { return "fermion"; } 00127 00128 }; 00129 00130 /** \brief Fermion with finite-temperature thermodynamics 00131 [abstract base] 00132 00133 This is an abstract base for the computation of 00134 finite-temperature fermionic statistics. Different children 00135 (e.g. \ref eff_fermion and \ref rel_fermion) use different 00136 techniques to computing the momentum integrations. 00137 00138 Because massless fermions at finite temperature are much 00139 simpler, there are separate member functions included in this 00140 class to handle them. The functions massless_calc_density() and 00141 massless_calc_mu() compute the thermodynamics of massless 00142 fermions at finite temperature given the density or the chemical 00143 potentials. The functions massless_pair_density() and 00144 massless_pair_mu() perform the same task, but automatically 00145 include antiparticles. 00146 00147 The function massless_calc_density() uses a \ref root object to 00148 solve for the chemical potential as a function of the density. 00149 The default is an object of type cern_mroot_root. The function 00150 massless_pair_density() does not need to use the \ref root 00151 object because of the simplification afforded by the inclusion 00152 of antiparticles. 00153 00154 \future Create a Chebyshev approximation for inverting the 00155 the Fermi functions for massless_calc_density() functions? 00156 00157 \htmlonly 00158 This Mathematica notebook contains the derivations of related 00159 series expansions and some algebra for the massless_pair() 00160 functions. 00161 <a href="../extras/fermion.nb"> 00162 fermion.nb</a>, and 00163 <a href="../extras/fermion.pdf"> 00164 fermion.pdf</a>. 00165 \endhtmlonly 00166 \latexonly 00167 This Mathematica notebook contains the derivations of related 00168 series expansions and some algebra for the massless\_pair() 00169 functions. 00170 \begin{verbatim} 00171 doc/o2scl/extras/fermion.nb 00172 doc/o2scl/extras/fermion.pdf 00173 \end{verbatim} 00174 \endlatexonly 00175 */ 00176 class fermion_T : public fermion { 00177 00178 public: 00179 00180 /// Create a fermion with mass \c mass and degeneracy \c dof. 00181 fermion_T(double mass=0, double dof=0); 00182 00183 virtual ~fermion_T() { 00184 } 00185 00186 /** 00187 \brief Calculate properties as function of chemical potential 00188 */ 00189 virtual int calc_mu(const double temper)=0; 00190 00191 /** \brief Calculate properties as function of density 00192 */ 00193 virtual int calc_density(const double temper)=0; 00194 00195 /** \brief Calculate properties with antiparticles as function of 00196 chemical potential 00197 */ 00198 virtual int pair_mu(const double temper)=0; 00199 00200 /** \brief Calculate properties with antiparticles as function of 00201 density 00202 */ 00203 virtual int pair_density(const double temper)=0; 00204 00205 #ifdef O2SCL_NEVER_DEFINED 00206 /** 00207 \brief Degenerate expansion for specific heat 00208 00209 This is a temporary location and is also unchecked. 00210 */ 00211 double deg_specific_heat(double T) { 00212 double ret; 00213 if (non_interacting) { 00214 nu=mu; 00215 ms=m; 00216 } 00217 double sqd=nu*nu-ms*ms; 00218 ret=o2scl_const::pi2*T*nu/sqd* 00219 (1.0-o2scl_const::pi2*T*T* 00220 (5.0*pow(ms,4.0)+4.0*ms*ms*nu*nu+14.0*pow(nu,4.0))/ 00221 15.0/nu/nu/sqd/sqd); 00222 return ret; 00223 } 00224 #endif 00225 00226 /// \name Massless fermions 00227 //@{ 00228 /// Finite temperature massless fermions 00229 virtual int massless_calc_mu(const double temper); 00230 00231 /// Finite temperature massless fermions 00232 virtual int massless_calc_density(const double temper); 00233 00234 /** 00235 \brief Finite temperature massless fermions and antifermions 00236 */ 00237 int massless_pair_mu(const double temper); 00238 00239 /** 00240 \brief Finite temperature massless fermions and antifermions 00241 00242 In the cases \f$ n^3 >> T \f$ and \f$ T >> n^3 \f$ , 00243 expansions are used instead of the exact formulas to avoid 00244 loss of precision. 00245 00246 In particular, using the parameter 00247 \f[ 00248 \alpha = \frac{g^2 \pi^2 T^6}{243 n^2} 00249 \f] 00250 and defining the expression 00251 \f[ 00252 \mathrm{cbt} = \alpha^{-1/6} \left( -1 + \sqrt{1+\alpha}\right)^{1/3} 00253 \f] 00254 we can write the chemical potential as 00255 \f[ 00256 \mu = \frac{\pi T}{\sqrt{3}} \left(\frac{1}{\mathrm{cbt}} - 00257 \mathrm{cbt} \right) 00258 \f] 00259 00260 These expressions, however, do not work well when \f$ \alpha 00261 \f$ is very large or very small, so series expansions are 00262 used whenever \f$ \alpha > 10^{4} \f$ or 00263 \f$ \alpha < 3 \times 10^{-4} \f$. For large \f$ \alpha \f$, 00264 \f[ 00265 \left(\frac{1}{\mathrm{cbt}} - 00266 \mathrm{cbt} \right) \approx 00267 \frac{2^{1/3}}{\alpha^{1/6}} - 00268 \frac{\alpha^{1/6}}{2^{1/3}} + 00269 \frac{\alpha^{5/6}}{6~2^{2/3}} + 00270 \frac{\alpha^{7/6}}{12~2^{1/3}} - 00271 \frac{\alpha^{11/6}}{18~2^{2/3}} - 00272 \frac{5 \alpha^{13/6}}{144~2^{1/3}} + 00273 \frac{77 \alpha^{17/6}}{2592~2^{2/3}} 00274 \f] 00275 and for small \f$ \alpha \f$, 00276 \f[ 00277 \left(\frac{1}{\mathrm{cbt}} - 00278 \mathrm{cbt} \right) \approx 00279 \frac{2}{3} \sqrt{\frac{1}{\alpha}} - 00280 \frac{8}{81} \left(\frac{1}{\alpha}\right)^{3/2} + 00281 \frac{32}{729} \left(\frac{1}{\alpha}\right)^{5/2} 00282 \f] 00283 00284 This approach works to within about 1 \part in \f$ 10^{12} \f$, 00285 and is tested in <tt>fermion_ts.cpp</tt>. 00286 00287 \future This could be improved by including more terms 00288 in the expansions. 00289 00290 */ 00291 int massless_pair_density(const double temper); 00292 //@} 00293 00294 /** \brief Set the solver for use in massless_calc_density() */ 00295 int set_massless_root(root<const double, funct<const double> > &rp) { 00296 massless_root=&rp; 00297 return 0; 00298 } 00299 00300 /** 00301 \brief The default solver for massless_calc_density() 00302 00303 We default to cern_mroot_root here since we don't have 00304 a bracket or a derivative. 00305 */ 00306 cern_mroot_root<const double, funct<const double> > def_massless_root; 00307 00308 /// Return string denoting type ("fermion_T") 00309 virtual const char *type() { return "fermion_T"; } 00310 00311 #ifndef DOXYGENP 00312 00313 protected: 00314 00315 /// A pointer to the solver for massless fermions 00316 root<const double, funct<const double> > *massless_root; 00317 00318 /// The function to solve for massless fermions 00319 int massless_fun(double x, double &y, const double &temper); 00320 00321 #endif 00322 00323 }; 00324 00325 template<> int io_tlate<fermion>::input 00326 (cinput *co, in_file_format *ins, fermion *f); 00327 template<> int io_tlate<fermion>::output 00328 (coutput *co, out_file_format *outs, fermion *f); 00329 template<> const char *io_tlate<fermion>::type(); 00330 00331 typedef io_tlate<fermion> fermion_io_type; 00332 00333 #ifndef DOXYGENP 00334 } 00335 #endif 00336 00337 #endif
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page