#include <cfl_njl_eos.h>
The original Lagrangian is
where is the quark number chemical potential. couplings
,
, and
ultra-violet three-momentum cutoff,
The thermodynamic potential is
where runs over all nine (three colors times three flavors) quarks. We assume that the condensates are independent of color and that the quark chemical potentials are of the form
with
With these assumptions, the thermodynamic potential as given by the function thd_potential(), is a function of 12 variables
The individual terms are
where are the eigenvalues of the (72 by 72) matrix (calculated by the function eigenvalues())
and is the charge conjugation matrix (in the Dirac representation).
The values of the various condensates are usually determined by the condition
Note that setting fixed_mass to true
and setting all of the gaps to zero when gap_limit
is less than zero will reproduce an analog of the bag model with a momentum cutoff.
The variable nambujl_eos::fromqq is automatically set to true in the constructor, as computations with fromqq=false
are not implemented.
Created for Steiner02.
Definition at line 208 of file cfl_njl_eos.h.
Public Member Functions | |
virtual int | set_parameters (double lambda=0.0, double fourferm=0.0, double sixferm=0.0, double fourgap=0.0) |
Set the parameters and the bag constant 'B0'. | |
virtual int | calc_eq_temp_p (quark &u, quark &d, quark &s, double &qq1, double &qq2, double &qq3, double &gap1, double &gap2, double &gap3, double mu3, double mu8, double &n3, double &n8, thermo &qb, double temper) |
Calculate the EOS. | |
virtual int | test_derivatives (double lmom, double mu3, double mu8, test_mgr &t) |
Check the derivatives specified by eigenvalues(). | |
virtual int | eigenvalues (double lmom, double mu3, double mu8, double egv[36], double dedmuu[36], double dedmud[36], double dedmus[36], double dedmu[36], double dedmd[36], double dedms[36], double dedu[36], double dedd[36], double deds[36], double dedmu3[36], double dedmu8[36]) |
Calculate the energy eigenvalues as a function of the momentum. | |
int | set_quartic (quartic_real_coeff &q) |
Set the routine for solving quartics. | |
int | test_integration (test_mgr &t) |
Test the integration routines. | |
int | test_normal_eigenvalues (test_mgr &t) |
Test the routine to compute the eigenvalues of non-superfluid fermions. | |
int | test_gapped_eigenvalues (test_mgr &t) |
Test the routine to compute the eigenvalues of superfluid fermions. | |
virtual const char * | type () |
Return string denoting type ("cfl_njl_eos"). | |
Data Fields | |
double | eq_limit |
The equal mass threshold. | |
bool | integ_test |
Set to true to test the integration (default false). | |
cern_quartic_real_coeff | def_quartic |
The default quartic routine. | |
double | gap_limit |
Smallest allowable gap (default 0.0). | |
bool | zerot |
If this is true, then finite temperature corrections are ignored (default false). | |
bool | fixed_mass |
Use a fixed quark mass and ignore the quark condensates. | |
bool | color_neut |
If true, then ensure color neutrality. | |
double | GD |
Diquark coupling constant (default 3 G/4). | |
double | inte_epsabs |
The absolute precision for the integration (default ![]() | |
double | inte_epsrel |
The relative precision for the integration (default ![]() | |
size_t | inte_npoints |
The number of points used in the last integration (default 0). | |
Protected Member Functions | |
virtual int | integrands (double p, double res[]) |
The integrands. | |
int | normal_eigenvalues (double m, double lmom, double mu, double lam[2], double dldmu[2], double dldm[2]) |
Compute ungapped eigenvalues and the appropriate derivatives. | |
int | gapped_eigenvalues (double m1, double m2, double lmom, double mu1, double mu2, double tdelta, double lam[4], double dldmu1[4], double dldmu2[4], double dldm1[4], double dldm2[4], double dldg[4]) |
Treat the simply gapped quarks in all cases gracefully. | |
For the integration | |
double | rescale_error (double err, double result_abs, double result_asc) |
The error scaling function for integ_err. | |
int | integ_err (double a, double b, const size_t nr, ovector &res, double &err2) |
A new version of gsl_inte_qng to integrate several functions at the same time. | |
Protected Attributes | |
double | temper |
Temperature. | |
double | smu3 |
3rd gluon chemical potential | |
double | smu8 |
8th gluon chemical potential | |
Numerical methods | |
quartic_real_coeff * | quartic |
The routine to solve quartics. | |
For computing eigenvalues | |
omatrix_cx | iprop |
Inverse propagator matrix. | |
omatrix_cx | eivec |
The eigenvectors. | |
omatrix_cx | dipdgapu |
The derivative of the inverse propagator wrt the ds gap. | |
omatrix_cx | dipdgapd |
The derivative of the inverse propagator wrt the us gap. | |
omatrix_cx | dipdgaps |
The derivative of the inverse propagator wrt the ud gap. | |
ovector | eval |
The eigenvalues. | |
gsl_eigen_hermv_workspace * | w |
Workspace for eigenvalue computation. | |
Private Member Functions | |
cfl_njl_eos (const cfl_njl_eos &) | |
cfl_njl_eos & | operator= (const cfl_njl_eos &) |
virtual int calc_eq_temp_p | ( | quark & | u, | |
quark & | d, | |||
quark & | s, | |||
double & | qq1, | |||
double & | qq2, | |||
double & | qq3, | |||
double & | gap1, | |||
double & | gap2, | |||
double & | gap3, | |||
double | mu3, | |||
double | mu8, | |||
double & | n3, | |||
double & | n8, | |||
thermo & | qb, | |||
double | temper | |||
) | [virtual] |
Calculate the EOS from the quark condensates in u.qq
, d.qq
and s.qq
. Return the mass gap equations in qq1
, qq2
, qq3
, and the normal gap equations in gap1
, gap2
, and gap3
.
Using fromqq=false
as in nambujl_eos and nambujl_eos does not work here and will return an error. Also, the quarks must be set through quark_eos::quark_set() before use.
If all of the gaps are less than gap_limit, then the nambujl_eos::calc_temp_p() is used, and gap1
, gap2
, and gap3
are set to equal u.del
, d.del
, and s.del
, respectively.
Reimplemented in cfl6_eos.
virtual int eigenvalues | ( | double | lmom, | |
double | mu3, | |||
double | mu8, | |||
double | egv[36], | |||
double | dedmuu[36], | |||
double | dedmud[36], | |||
double | dedmus[36], | |||
double | dedmu[36], | |||
double | dedmd[36], | |||
double | dedms[36], | |||
double | dedu[36], | |||
double | dedd[36], | |||
double | deds[36], | |||
double | dedmu3[36], | |||
double | dedmu8[36] | |||
) | [virtual] |
Given the momentum mom
, and the chemical potentials associated with the third and eighth gluons (mu3
and mu8
), the energy eigenvalues are computed in egv[0] ... egv[35].
int gapped_eigenvalues | ( | double | m1, | |
double | m2, | |||
double | lmom, | |||
double | mu1, | |||
double | mu2, | |||
double | tdelta, | |||
double | lam[4], | |||
double | dldmu1[4], | |||
double | dldmu2[4], | |||
double | dldm1[4], | |||
double | dldm2[4], | |||
double | dldg[4] | |||
) | [protected] |
This function uses the quarks q1
and q2
to construct the eigenvalues of the inverse propagator, properly handling the either zero or finite quark mass and either zero or finite quark gaps. In the case of finite quark mass and finite quark gaps, the quartic solver is used.
The chemical potentials are separated so we can add the color chemical potentials to the quark chemical potentials if necessary.
This function is used by eigenvalues(). It does not work for the "ur-dg-sb" set of quarks which are paired in a non-trivial way.
virtual int integrands | ( | double | p, | |
double | res[] | |||
) | [protected, virtual] |
Reimplemented in cfl6_eos.
virtual int set_parameters | ( | double | lambda = 0.0 , |
|
double | fourferm = 0.0 , |
|||
double | sixferm = 0.0 , |
|||
double | fourgap = 0.0 | |||
) | [virtual] |
This function allows the user to specify the momentum cutoff, lambda
, the four-fermion coupling fourferm
, the six-fermion coupling from the 't Hooft interaction sixferm
, and the color-superconducting coupling, fourgap
. If 0.0 is given for any of the values, then the default is used ().
If the four-fermion coupling that produces a gap is not specified, it is automatically set to 3/4 G, which is the value obtained from the Fierz transformation.
The value of the shift in the bag constant nambujl_eos::B0 is automatically calculated to ensure that the vacuum has zero energy density and zero pressure. The functions set_quarks() and set_thermo() must be used before hand to specify the quark and thermo objects.
Slightly better accuracy (with slower execution times) can be achieved using gsl_poly_real_coeff which polishes the roots of the quartics. For example
cfl_njl_eos cfl; gsl_poly_real_coeff gp; cfl.set_quartic(gp);
Definition at line 318 of file cfl_njl_eos.h.
double gap_limit |
If any of the gaps are below this value, then it is assumed that they are zero and the equation of state is simplified accordingly. If all of the gaps are less than gap_limit, then the results from nambujl_eos are used in calc_eq_temp_p(), calc_temp_p() and thd_potential().
Definition at line 342 of file cfl_njl_eos.h.
double GD |
The default value is the one derived from a Fierz transformation. (Buballa04)
Definition at line 368 of file cfl_njl_eos.h.
double inte_epsabs |
This is analogous to gsl_inte::epsabs
Definition at line 376 of file cfl_njl_eos.h.
double inte_epsrel |
This is analogous to gsl_inte::epsrel
Definition at line 384 of file cfl_njl_eos.h.
size_t inte_npoints |
This returns 21, 43, or 87 depending on the number of function evaluations needed to obtain the desired precision. If it the routine failes to obtain the desired precision, then this variable is set to 88.
Definition at line 395 of file cfl_njl_eos.h.
bool zerot |
This implements some simplifications in the momentum integration that are not possible at finite temperature.
Definition at line 351 of file cfl_njl_eos.h.
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