#include <cfl_njl_eos_old.h>
The variable B0 must be set before use.
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.
Definition at line 202 of file cfl_njl_eos_old.h.
Numerical methods | |
gsl_mmin_conf< void *, multi_funct< void * > > | def_min1 |
The first default minimizer. | |
gsl_mmin_conf< void *, multi_funct< void * > > | def_min2 |
The second default minimizers. | |
cern_quartic_real_coeff | def_quartic |
The default quartic routine. | |
int | set_minimizers (multi_min< void *, multi_funct< void * > > &mi1, multi_min< void *, multi_funct< void * > > &mi2) |
Set the minimizers for the thermodynamic potential. | |
int | set_quartic (quartic_real_coeff &q) |
Set the routine for solving quartics. | |
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_temp_p (quark &u, quark &d, quark &s, const double temper, thermo &th) |
Calculate equation of state as a function of chemical potentials. | |
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, const double temper) |
Calculate the EOS. | |
virtual int | eigenvalues (double mom, ovector_view &egv, double mu3, double mu8) |
Calculate the energy eigenvalues as a function of the momentum. | |
virtual double | thd_potential (quark &u, quark &d, quark &s, double mu3, double mu8, const double ltemper) |
Direct calculation of the thermodynamic potential. | |
virtual const char * | type () |
Return string denoting type ("cfl_njl_eos_old"). | |
virtual double | tpot (double var, void *&pa) |
The function used to take derivatives of the thermodynamic potential (used by calc_temp_p() ). | |
Data Fields | |
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). | |
bool | faster |
Faster eigenvalues. | |
double | solver_qual |
The quality of the most recently obtained solution. | |
bool | use_solver |
If true, use the solver, otherwise, use the minimizers. | |
Protected Member Functions | |
int | solve_fun_p (size_t nv, const ovector_view &x, ovector_view &y, void *&vp) |
Function specifying gap equations to solve. | |
int | min_fun2 (size_t nv, const ovector_view &x, double &y, void *&vp) |
Outer minimization function for color neutrality. | |
int | min_fun (size_t nv, const ovector_view &x, double &y, void *&vp) |
Inner minimization function for mass and gap equations. | |
double | derivfn (double *dvar) |
A function that takes the derivative of the thermodynamic potential with respect to the variable whose address is 'dvar'. | |
double | tpot_integrand (double p, void *&pa) |
The integrand of the thermodynamic potential. | |
double | entropy_integrand (double p, void *&pa) |
The integrand of the entropy. | |
int | gapped_eigenvalues (quark *q1, quark *q2, double mom, double mu1, double mu2, double delta, ovector_cx_view &ev) |
Treat the simply gapped quarks in all cases gracefully. | |
Protected Attributes | |
double | temper |
Temperature. | |
double | smu3 |
3rd gluon chemical potential | |
double | smu8 |
8th gluon chemical potential | |
double * | var |
Pointer to the variable to be differentiated by. | |
Numerical methods | |
quartic_real_coeff * | quartic |
The routine to solve quartics. | |
bool | quartic_set |
true if the quartic routine has been set | |
gsl_deriv< void *, funct< void * > > | de |
The derivative routine. | |
multi_min< void *, multi_funct < void * > > * | min1 |
Pointer to inner minimizer. | |
multi_min< void *, multi_funct < void * > > * | min2 |
Pointer to outer minimizer. | |
For computing eigenvalues | |
omatrix_cx | m |
Inverse propagator matrix. | |
omatrix_cx | m2 |
Inverse propagator matrix. | |
ovector | eval |
Storage for eigenvalues. | |
ovector | eval2 |
Storage for eigenvalues. | |
gsl_eigen_herm_workspace * | w |
Workspace for eigenvalue computation. | |
gsl_eigen_herm_workspace * | w2 |
Workspace for eigenvalue computation. | |
bool | ugap_solve |
Switches for solving the gap equations. | |
bool | dgap_solve |
bool | sgap_solve |
Private Member Functions | |
cfl_njl_eos_old (const cfl_njl_eos_old &) | |
cfl_njl_eos_old & | operator= (const cfl_njl_eos_old &) |
virtual int set_parameters | ( | double | lambda = 0.0 , |
|
double | fourferm = 0.0 , |
|||
double | sixferm = 0.0 , |
|||
double | fourgap = 0.0 | |||
) | [virtual] |
Set the parameters and the bag constant 'B0'.
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.
virtual int calc_temp_p | ( | quark & | u, | |
quark & | d, | |||
quark & | s, | |||
const double | temper, | |||
thermo & | th | |||
) | [virtual] |
Calculate equation of state as a function of chemical potentials.
If any of the gaps given in quark::del are greater than gap_limit, then the respective gap equations are solved. If fixed_mass is false, then the quark condensates are solved for. If color_neut is true, then color neutrality is enforced.
Reimplemented from nambujl_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, | |||
const double | temper | |||
) | [virtual] |
Calculate the EOS.
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_old.
virtual int eigenvalues | ( | double | mom, | |
ovector_view & | egv, | |||
double | mu3, | |||
double | mu8 | |||
) | [virtual] |
Calculate the energy eigenvalues as a function of the momentum.
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]. No space is allocated for the array by the function.
Reimplemented in cfl6_eos_old.
virtual double tpot | ( | double | var, | |
void *& | pa | |||
) | [virtual] |
The function used to take derivatives of the thermodynamic potential (used by calc_temp_p() ).
To use this function just to compute the thermodynamic potential, make sure var is 0 (or NULL
) and use
tpot(0.0,0);
Reimplemented in cfl6_eos_old.
int gapped_eigenvalues | ( | quark * | q1, | |
quark * | q2, | |||
double | mom, | |||
double | mu1, | |||
double | mu2, | |||
double | delta, | |||
ovector_cx_view & | ev | |||
) | [protected] |
Treat the simply gapped quarks in all cases gracefully.
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.
The default quartic routine.
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_old cfl; gsl_poly_real_coeff gp; cfl.set_quartic(gp);
Definition at line 324 of file cfl_njl_eos_old.h.
double gap_limit |
Smallest allowable gap (default 0.0).
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 336 of file cfl_njl_eos_old.h.
bool zerot |
If this is true, then finite temperature corrections are ignored (default false).
This implements some simplifications in the momentum integration that are not possible at finite temperature.
Definition at line 344 of file cfl_njl_eos_old.h.
double GD |
Diquark coupling constant (default 3 G/4).
The default value is the one derived from a Fierz transformation. (Buballa04)
Definition at line 360 of file cfl_njl_eos_old.h.
bool faster |
Faster eigenvalues.
This assumes that the strange quark mass can be used as an effective chemical potential and calculates the eigenvalues to first order in the four variables
Definition at line 384 of file cfl_njl_eos_old.h.
bool quartic_set [protected] |
true if the quartic routine has been set
This is used for the I/O.
Definition at line 480 of file cfl_njl_eos_old.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