23 #ifndef O2SCL_FERMION_DERIV_REL_H
24 #define O2SCL_FERMION_DERIV_REL_H
34 #include <o2scl/constants.h>
35 #include <o2scl/root_cern.h>
36 #include <o2scl/inte.h>
37 #include <o2scl/inte_qag_gsl.h>
38 #include <o2scl/inte_qagiu_gsl.h>
40 #include <o2scl/part_deriv.h>
41 #include <o2scl/fermion_rel.h>
43 #ifndef DOXYGEN_NO_O2NS
250 template<
class fp_t=
double>
362 fp_t prefac=f.
g/2.0/this->
pi2;
417 funct density_T_fun_f=
420 this,std::placeholders::_1,std::ref(f),temper);
423 O2SCL_ERR2(
"dndT integration (ndeg) failed in ",
424 "fermion_deriv_rel::calc_mu().",
430 funct density_mu_fun_f=
433 this,std::placeholders::_1,std::ref(f),temper);
436 O2SCL_ERR2(
"dndmu integration (ndeg) failed in ",
437 "fermion_deriv_rel::calc_mu().",
443 funct entropy_T_fun_f=
446 this,std::placeholders::_1,std::ref(f),temper);
449 O2SCL_ERR2(
"dsdT integration (ndeg) failed in ",
450 "fermion_deriv_rel_tl<fp_t>::calc_mu().",
exc_efailed);
470 O2SCL_ERR2(
"Zero density in degenerate limit in fermion_deriv_rel::",
471 "calc_mu(). Variable deg_limit set improperly?",
510 funct deg_density_mu_fun_f=
513 this,std::placeholders::_1,std::ref(f),temper);
522 O2SCL_ERR2(
"dndmu integration (deg) failed in ",
523 "fermion_deriv_rel_tl<fp_t>::calc_mu().",
exc_efailed);
528 funct deg_density_T_fun_f=std::bind
531 this,std::placeholders::_1,std::ref(f),temper);
538 O2SCL_ERR2(
"dndT integration (deg) failed in ",
539 "fermion_deriv_rel_tl<fp_t>::calc_mu().",
exc_efailed);
544 funct deg_entropy_T_fun_f=std::bind
547 this,std::placeholders::_1,std::ref(f),temper);
554 O2SCL_ERR2(
"dsdT integration (deg) failed in ",
555 "fermion_deriv_rel_tl<fp_t>::calc_mu().",
exc_efailed);
562 if (!std::isfinite(f.
en)) {
564 "fermion_deriv_rel_tl<fp_t>::calc_mu().",
exc_efailed);
666 virtual const char *
type() {
return "fermion_deriv_rel"; };
670 #ifndef DOXYGEN_NO_O2NS
689 fp_t k=u*(T), E, ret;
694 ret=k*k*(E-f.
nu)/T*ff*(1.0-ff);
696 ret=(2.0*k*k/T+E*E/T-E*(f.
nu)/T-k*k*(f.
nu)/T/E)*
704 ret=k*k*(E-f.
nu)/T*ff*(1.0-ff);
706 ret=(2.0*k*k/T+E*E/T-E*(f.
nu+f.
m)/T-k*k*(f.
nu+f.
m)/T/E)*
714 fp_t k=u*(T), E, ret;
736 fp_t entropy_T_fun(fp_t u, fermion_deriv &f, fp_t T) {
738 if (f.inc_rest_mass) {
742 ret=T*k*k*ff*(1.0-ff)*pow(E-f.nu,2.0)/pow(T,3.0);
745 (pow(E,3.0)+3.0*E*k*k-(E*E+k*k)*(f.nu))*
753 ret=T*k*k*ff*(1.0-ff)*pow(E-f.nu,2.0)/pow(T,3.0);
755 ret=(E-f.m-f.nu)/E/T*
756 (pow(E,3.0)+3.0*E*k*k-(E*E+k*k)*(f.nu+f.m))*
763 fp_t density_ms_fun(fp_t u, fermion_deriv &f, fp_t T) {
765 if (f.inc_rest_mass) {
769 ret=-k*k*f.ms/(E)/T*ff*(1.0-ff);
778 ret=-k*k*f.ms/(E+f.m)/T*ff*(1.0-ff);
793 fp_t deg_density_T_fun(fp_t k, fermion_deriv &f, fp_t T) {
795 if (f.inc_rest_mass) {
799 ret=k*k*(E-f.nu)/T/T*ff*(1.0-ff);
801 ret=(2.0*k*k/T+E*E/T-E*(f.nu)/T-k*k*(f.nu)/T/E)*
809 ret=k*k*(E-f.nu)/T/T*ff*(1.0-ff);
811 ret=(2.0*k*k/T+E*E/T-E*(f.nu+f.m)/T-k*k*(f.nu+f.m)/T/E)*
818 fp_t deg_density_mu_fun(fp_t k, fermion_deriv &f, fp_t T) {
820 if (f.inc_rest_mass) {
824 ret=k*k/T*ff*(1.0-ff);
833 ret=k*k/T*ff*(1.0-ff);
841 fp_t deg_entropy_T_fun(fp_t k, fermion_deriv &f, fp_t T) {
844 if (f.inc_rest_mass) {
847 ret=k*k*ff*(1.0-ff)*pow(E-f.nu,2.0)/pow(T,3.0);
850 (pow(E,3.0)+3.0*E*k*k-(E*E+k*k)*f.nu)*ff;
855 ret=k*k*ff*(1.0-ff)*pow(E-f.nu-f.m,2.0)/pow(T,3.0);
857 ret=(E-f.m-f.nu)/E/T/T*
858 (pow(E,3.0)+3.0*E*k*k-(E*E+k*k)*(f.nu+f.m))*ff;
864 fp_t deg_density_ms_fun(fp_t k, fermion_deriv &f, fp_t T) {
866 if (f.inc_rest_mass) {
870 ret=-k*k*f.ms/E/T*ff*(1.0-ff);
879 ret=-k*k*f.ms/(E+f.m)/T*ff*(1.0-ff);
898 #ifndef DOXYGEN_NO_O2NS