Fermion with finite-temperature thermodynamics [abstract base]. More...
#include <fermion.h>
This is an abstract base for the computation of finite-temperature fermionic statistics. Different children (e.g. fermion_eff and fermion_rel) use different techniques to computing the momentum integrations.
Because massless fermions at finite temperature are much simpler, there are separate member functions included in this class to handle them. The functions massless_calc_density() and massless_calc_mu() compute the thermodynamics of massless fermions at finite temperature given the density or the chemical potentials. The functions massless_pair_density() and massless_pair_mu() perform the same task, but automatically include antiparticles.
The function massless_calc_density() uses a root object to solve for the chemical potential as a function of the density. The default is an object of type root_cern. The function massless_pair_density() does not need to use the root object because of the simplification afforded by the inclusion of antiparticles.
Public Member Functions | |
virtual bool | calc_mu_ndeg (fermion &f, double temper, double prec=1.0e-18, bool inc_antip=false) |
Non-degenerate expansion for fermions. More... | |
virtual bool | calc_mu_deg (fermion &f, double temper, double prec=1.0e-18) |
Degenerate expansion for fermions. More... | |
virtual void | calc_mu (fermion &f, double temper)=0 |
Calculate properties as function of chemical potential. | |
virtual int | calc_density (fermion &f, double temper)=0 |
Calculate properties as function of density. More... | |
virtual void | pair_mu (fermion &f, double temper)=0 |
Calculate properties with antiparticles as function of chemical potential. | |
virtual int | pair_density (fermion &f, double temper)=0 |
Calculate properties with antiparticles as function of density. More... | |
void | set_massless_root (root<> &rp) |
Set the solver for use in massless_calc_density() | |
virtual const char * | type () |
Return string denoting type ("fermion_thermo") | |
template<class fermion_t > | |
bool | calc_mu_deg_tlate (fermion_t &f, double temper, double prec) |
Calculate thermodynamic properties from the chemical potential using a degenerate expansion. | |
void | ndeg_terms (size_t j, double tt, double xx, double m, bool inc_rest_mass, bool inc_antip, double &pterm, double &nterm, double &enterm) |
Compute a term in the nondegenerate expansion. | |
template<class fermion_t > | |
bool | calc_mu_ndeg_tlate (fermion_t &f, double temper, double prec, bool inc_antip) |
Calculate thermodynamic properties from the chemical potential using a nondegenerate expansion. | |
Massless fermions | |
virtual void | massless_calc_mu (fermion &f, double temper) |
Finite temperature massless fermions. | |
virtual void | massless_calc_density (fermion &f, double temper) |
Finite temperature massless fermions. | |
virtual void | massless_pair_mu (fermion &f, double temper) |
Finite temperature massless fermions and antifermions. | |
virtual void | massless_pair_density (fermion &f, double temper) |
Finite temperature massless fermions and antifermions. More... | |
![]() | |
void | kf_from_density (fermion &f) |
Calculate the Fermi momentum from the density. More... | |
void | energy_density_zerot (fermion &f) |
Energy density at T=0 from fermion::kf and part::ms. More... | |
void | pressure_zerot (fermion &f) |
Pressure at T=0 from fermion::kf and part::ms. More... | |
virtual void | calc_mu_zerot (fermion &f) |
Zero temperature fermions from part::mu or part::nu and part::ms. | |
virtual void | calc_density_zerot (fermion &f) |
Zero temperature fermions from part::n and part::ms. | |
Public Attributes | |
root_cern | def_massless_root |
The default solver for massless_calc_density() More... | |
Protected Member Functions | |
double | massless_solve_fun (double x, fermion &f, double temper) |
Solve for the chemical potential for massless fermions. | |
Protected Attributes | |
root * | massless_root |
A pointer to the solver for massless fermions. | |
|
pure virtual |
Implemented in o2scl::fermion_rel, o2scl::fermion_eff, and o2scl::fermion_nonrel.
|
virtual |
Attempts to evaulate thermodynamics of a degenerate fermion. If the result is accurate to within the requested precision, this function returns true
, and otherwise this function returns false
and the values in stored in the pr
, n
, en
, and ed
field are meaningless.
The pressure, density, and energy density, should be accurate to the requested precision, but the first term in the series expansion for the entropy is zero, so the entropy is one order lower in accuracy.
|
virtual |
Attempts to evaluate thermodynamics of a non-degenerate fermion. If the result is accurate to within the requested precision, this function returns true
, and otherwise this function returns false
and the values in stored in the pr
, n
, en
, and ed
field are meaningless.
If is negative and sufficiently far from zero, then the thermodynamic quantities are smaller than the smallest representable double-precision number. In this case, this function will return
true
and report all quantities as zero.
Defining ,
, and
the pressure in the non-degenerate limit (
) is (Johns96)
where
The density is then
and the entropy density is
This function is accurate over a wide range of conditions when .
The ratio of the nth term to the first term in the pressure series is
This function currently uses 20 terms in the series and immediately returns false
if is greater than
prec
In the nondegenerate and nonrelativistic ( ) limit, the argument to the Bessel functions and the exponential becomes too large. In this case, it's better to use the expansions, e.g. for
,
The current code currently goes up to in the expansion, which is enough for the default precision of
since
.
|
virtual |
In the cases and
, expansions are used instead of the exact formulas to avoid loss of precision.
In particular, using the parameter
and defining the expression
we can write the chemical potential as
These expressions, however, do not work well when is very large or very small, so series expansions are used whenever
or
. For small
,
and for large ,
This approach works to within about 1 part in , and is tested in
fermion_ts.cpp
.
|
pure virtual |
Implemented in o2scl::fermion_rel, o2scl::fermion_eff, and o2scl::fermion_nonrel.
root_cern o2scl::fermion_thermo::def_massless_root |
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).