Go to the documentation of this file.
32 #include <boost/numeric/ublas/vector.hpp>
33 #include <boost/numeric/ublas/vector_proxy.hpp>
34 #include <boost/numeric/ublas/matrix.hpp>
35 #include <boost/numeric/ublas/matrix_proxy.hpp>
37 #include <gsl/gsl_math.h>
38 #include <gsl/gsl_eigen.h>
39 #include <gsl/gsl_complex_math.h>
40 #include <gsl/gsl_complex.h>
41 #include <gsl/gsl_sort_vector.h>
42 #include <gsl/gsl_poly.h>
44 #include <o2scl/constants.h>
45 #include <o2scl/err_hnd.h>
46 #include <o2scl/mm_funct.h>
47 #include <o2scl/inte_qng_gsl.h>
48 #include <o2scl/poly.h>
49 #include <o2scl/test_mgr.h>
50 #include <o2scl/columnify.h>
52 #include <o2scl/part.h>
53 #include <o2scl/eos_quark_njl.h>
55 #ifndef DOXYGEN_NO_O2NS
236 virtual int set_parameters(
double lambda=0.0,
double fourferm=0.0,
237 double sixferm=0.0,
double fourgap=0.0);
259 double &qq2,
double &qq3,
double &gap1,
260 double &gap2,
double &gap3,
double mu3,
261 double mu8,
double &n3,
double &n8,
276 virtual int eigenvalues(
double lmom,
double mu3,
double mu8,
277 double egv[36],
double dedmuu[36],
278 double dedmud[36],
double dedmus[36],
279 double dedmu[36],
double dedmd[36],
280 double dedms[36],
double dedu[36],
281 double dedd[36],
double deds[36],
282 double dedmu3[36],
double dedmu8[36]);
385 virtual const char *
type() {
return "eos_quark_cfl"; };
387 #ifndef DOXYGEN_INTERNAL
407 virtual int integrands(
double p,
double res[]);
411 double lam[2],
double dldmu[2],
434 double mu1,
double mu2,
double tdelta,
435 double lam[4],
double dldmu1[4],
436 double dldmu2[4],
double dldm1[4],
437 double dldm2[4],
double dldg[4]);
478 gsl_eigen_hermv_workspace *
w;
490 int integ_err(
double a,
double b,
const size_t nr,
503 #ifndef DOXYGEN_NO_O2NS
gsl_eigen_hermv_workspace * w
Workspace for eigenvalue computation.
double inte_epsrel
The relative precision for the integration (default )
int integ_err(double a, double b, const size_t nr, ubvector &res, double &err2)
A new version of inte_qng_gsl to integrate several functions at the same time.
ubmatrix_complex dipdgapd
The derivative of the inverse propagator wrt the us gap.
virtual const char * type()
Return string denoting type ("eos_quark_cfl")
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.
virtual int test_derivatives(double lmom, double mu3, double mu8, test_mgr &t)
Check the derivatives specified by eigenvalues()
double gap_limit
Smallest allowable gap (default 0.0)
gsl_vector * eval
The eigenvalues.
double eq_limit
The equal mass threshold.
double smu3
3rd gluon chemical potential
bool integ_test
Set to true to test the integration (default false)
size_t inte_npoints
The number of points used in the last integration (default 0)
double rescale_error(double err, double result_abs, double result_asc)
The error scaling function for integ_err.
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 int integrands(double p, double res[])
The integrands.
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.
bool color_neut
If true, then ensure color neutrality.
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.
double smu8
8th gluon chemical potential
double temper
Temperature.
int test_integration(test_mgr &t)
Test the integration routines.
int set_quartic(quartic_real_coeff &q)
Set the routine for solving quartics.
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.
double inte_epsabs
The absolute precision for the integration (default )
ubmatrix_complex dipdgapu
The derivative of the inverse propagator wrt the ds gap.
gsl_matrix_complex * iprop
Inverse propagator matrix.
quartic_real_coeff_cern def_quartic
The default quartic routine.
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'.
double GD
Diquark coupling constant (default 3 G/4)
Nambu Jona-Lasinio model with a schematic CFL di-quark interaction at finite temperature.
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.
ubmatrix_complex dipdgaps
The derivative of the inverse propagator wrt the ud gap.
gsl_matrix_complex * eivec
The eigenvectors.
quartic_real_coeff * quartic
The routine to solve quartics.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).