fermion_deriv_nr.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2020, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_FERMION_DERIV_NR_H
24 #define O2SCL_FERMION_DERIV_NR_H
25 
26 /** \file fermion_deriv_nr.h
27  \brief File defining \ref o2scl::fermion_deriv_nr_tl
28 */
29 
30 #include <string>
31 #include <iostream>
32 #include <fstream>
33 #include <cmath>
34 #include <o2scl/constants.h>
35 #include <o2scl/root_cern.h>
36 
37 #include <o2scl/part_deriv.h>
38 #include <o2scl/classical_deriv.h>
39 
40 #ifndef DOXYGEN_NO_O2NS
41 namespace o2scl {
42 #endif
43 
44  /** \brief Equation of state for a nonrelativistic fermion
45 
46  This does not include the rest mass energy in the chemical
47  potential or the rest mass energy density in the energy density
48  to alleviate numerical precision problems at low densities
49 
50  This implements an equation of state for a nonrelativistic fermion
51  using direct integration. After subtracting the rest mass from
52  the chemical potentials, the distribution function is
53  \f[
54  \left\{1+\exp\left[\left(\frac{k^2}
55  {2 m^{*}}-\nu\right)/T\right]\right\}^{-1}
56  \f]
57  where \f$ \nu \f$ is the effective chemical potential, \f$ m \f$ is
58  the rest mass, and \f$ m^{*} \f$ is the effective mass.
59  For later use, we define \f$ E^{*} = k^2/2/m^{*} \f$ .
60 
61  Uncertainties are given in \ref unc.
62 
63  \b Evaluation \b of \b the \b derivatives
64 
65  The relevant derivatives of the distribution function are
66  \f[
67  \frac{\partial f}{\partial T}=
68  f(1-f)\frac{E^{*}-\nu}{T^2}
69  \f]
70  \f[
71  \frac{\partial f}{\partial \nu}=
72  f(1-f)\frac{1}{T}
73  \f]
74  \f[
75  \frac{\partial f}{\partial k}=
76  -f(1-f)\frac{k}{m^{*} T}
77  \f]
78  \f[
79  \frac{\partial f}{\partial m^{*}}=
80  f(1-f)\frac{k^2}{2 m^{*2} T}
81  \f]
82 
83  We also need the derivative of the entropy integrand w.r.t. the
84  distribution function, which is quite simple
85  \f[
86  {\cal S}\equiv f \ln f +(1-f) \ln (1-f) \qquad
87  \frac{\partial {\cal S}}{\partial f} = \ln
88  \left(\frac{f}{1-f}\right) =
89  \left(\frac{\nu-E^{*}}{T}\right)
90  \f]
91  where the entropy density is
92  \f[
93  s = - \frac{g}{2 \pi^2} \int_0^{\infty} {\cal S} k^2 d k
94  \f]
95 
96  The derivatives can be integrated directly
97  or they may be converted to integrals
98  over the distribution function through an integration by parts
99  \f[
100  \int_a^b f(k) \frac{d g(k)}{dk} dk = \left.f(k) g(k)\right|_{k=a}^{k=b}
101  - \int_a^b g(k) \frac{d f(k)}{dk} dk
102  \f]
103  using the distribution function for \f$ f(k) \f$ and 0 and \f$
104  \infty \f$ as the limits, we have
105  \f[
106  \frac{g}{2 \pi^2} \int_0^{\infty} \frac{d g(k)}{dk} f dk =
107  \frac{g}{2 \pi^2} \int_0^{\infty} g(k) f (1-f) \frac{k}{E^{*} T} dk
108  \f]
109  as long as \f$ g(k) \f$ vanishes at \f$ k=0 \f$ .
110  Rewriting,
111  \f[
112  \frac{g}{2 \pi^2} \int_0^{\infty} h(k) f (1-f) dk =
113  \frac{g}{2 \pi^2} \int_0^{\infty} f \frac{T m^{*}}{k}
114  \left[ h^{\prime} - \frac{h}{k}\right] d k
115  \f]
116  as long as \f$ h(k)/k \f$ vanishes at \f$ k=0 \f$ .
117 
118  \b Explicit \b forms
119 
120  1) The derivative of the density wrt the chemical potential
121  \f[
122  \left(\frac{d n}{d \mu}\right)_T =
123  \frac{g}{2 \pi^2} \int_0^{\infty} \frac{k^2}{T} f (1-f) dk
124  \f]
125  Using \f$ h(k)=k^2/T \f$ we get
126  \f[
127  \left(\frac{d n}{d \mu}\right)_T =
128  \frac{g}{2 \pi^2} \int_0^{\infty}
129  m^{*} f dk
130  \f]
131 
132  2) The derivative of the density wrt the temperature
133  \f[
134  \left(\frac{d n}{d T}\right)_{\mu} =
135  \frac{g}{2 \pi^2} \int_0^{\infty} \frac{k^2(E^{*}-\nu)}{T^2}
136  f (1-f) dk
137  \f]
138  Using \f$ h(k)=k^2(E^{*}-\nu)/T^2 \f$ we get
139  \f[
140  \left(\frac{d n}{d T}\right)_{\mu} =
141  \frac{g}{2 \pi^2} \int_0^{\infty} \frac{f}{T}
142  \left[m^{*} \left(E^{*}-\nu\right) -k^2\right] d k
143  \f]
144 
145  3) The derivative of the entropy wrt the chemical potential
146  \f[
147  \left(\frac{d s}{d \mu}\right)_T =
148  \frac{g}{2 \pi^2} \int_0^{\infty} k^2 f (1-f)
149  \frac{(E^{*}-\nu)}{T^2} dk
150  \f]
151  This verifies the Maxwell relation
152  \f[
153  \left(\frac{d s}{d \mu}\right)_T =
154  \left(\frac{d n}{d T}\right)_{\mu}
155  \f]
156 
157  4) The derivative of the entropy wrt the temperature
158  \f[
159  \left(\frac{d s}{d T}\right)_{\mu} =
160  \frac{g}{2 \pi^2} \int_0^{\infty} k^2 f (1-f)
161  \frac{(E^{*}-\nu)^2}{T^3} dk
162  \f]
163  Using \f$ h(k)=k^2 (E^{*}-\nu)^2/T^3 \f$
164  \f[
165  \left(\frac{d s}{d T}\right)_{\mu} =
166  \frac{g}{2 \pi^2} \int_0^{\infty}
167  f \frac{m^{*}}{T^2} \left[\left( E^{*}-\nu \right)^2
168  +\frac{2 k^2}{m^{*}} \left(E^{*}-\nu\right)\right] d k
169  \f]
170 
171  5) The derivative of the density wrt the effective mass
172  \f[
173  \left(\frac{d n}{d m^{*}}\right)_{T,\mu} =
174  \frac{g}{2 \pi^2} \int_0^{\infty}
175  \frac{k^2}{2 m^{* 2} T} f (1-f) k^2 dk
176  \f]
177  Using \f$ h(k)=k^4/(2 m^{* 2} T) \f$ we get
178  \f[
179  \left(\frac{d n}{d m^{*}}\right)_{T,\mu} =
180  \frac{g}{2 \pi^2} \int_0^{\infty} f
181  \frac{3 k^2}{2 m^{*}} d k
182  \f]
183 
184  <b>Conversion to unitless variables:</b>
185 
186  After integrating by parts
187  \f$ u = k^2/2/m^{*}/T \f$ and \f$ y=\mu/T \f$, so
188  \f[
189  k d k = m^{*} T d u
190  \f]
191  or
192  \f[
193  d k = \frac{m^{*} T}{\sqrt{2 m^{*} T u}} d u =
194  \sqrt{\frac{m^{*} T}{2 u}} d u
195  \f]
196 
197  1) The derivative of the density wrt the chemical potential
198  \f[
199  \left(\frac{d n}{d \mu}\right)_T =
200  \frac{g m^{* 3/2} \sqrt{T}}{2^{3/2} \pi^2} \int_0^{\infty}
201  u^{-1/2} f d u
202  \f]
203 
204  2) The derivative of the density wrt the temperature
205  \f[
206  \left(\frac{d n}{d T}\right)_{\mu} =
207  \frac{g m^{* 3/2} \sqrt{T}}
208  {2^{3/2} \pi^2} \int_0^{\infty} f d u
209  \left[ 3 u^{1/2} - y u^{-1/2}\right]
210  \f]
211 
212  4) The derivative of the entropy wrt the temperature
213  \f[
214  \left(\frac{d s}{d T}\right)_{\mu} =
215  \frac{g m^{* 3/2} T^{1/2}}{2^{3/2} \pi^2} \int_0^{\infty}
216  f \left[ 5 u^{3/2} - 6 y u^{1/2} + y^2 u^{-1/2}\right] d u
217  \f]
218 
219  5) The derivative of the density wrt the effective mass
220  \f[
221  \left(\frac{d n}{d m^{*}}\right)_{T,\mu} =
222  \frac{3 g m{* 1/2} T^{3/2}}{2^{3/2} \pi^2}
223  \int_0^{\infty} u^{1/2} f d u
224  \f]
225  */
226  template<class fp_t=double>
228 
229  public:
230 
231  /// Create a fermion with mass \c m and degeneracy \c g
233 
234  flimit=20.0;
235 
237  }
238 
239  virtual ~fermion_deriv_nr_tl() {
240  }
241 
242  /** \brief The limit for the Fermi functions (default 20.0)
243 
244  fermion_deriv_nr will ignore corrections smaller than about
245  \f$ \exp(-\mathrm{f{l}imit}) \f$ .
246  */
247  fp_t flimit;
248 
249  /// Storage for the most recently calculated uncertainties
251 
252  /** \brief Calculate properties as function of density
253  at \f$ T=0 \f$
254  */
256  if (f.non_interacting) { f.ms=f.m; }
257  f.kf=cbrt(6.0*this->pi2/f.g*f.n);
258  f.nu=f.kf*f.kf/2.0/f.ms;
259  f.ed=f.g*pow(f.kf,5.0)/20.0/this->pi2/f.ms;
260  if (f.inc_rest_mass) {
261  f.ed+=f.n*f.m;
262  f.nu+=f.m;
263  }
264  f.pr=-f.ed+f.n*f.nu;
265  f.en=0.0;
266 
267  if (f.non_interacting) { f.mu=f.nu; }
268 
269  f.dndT=0.0;
270  if (f.inc_rest_mass) {
271  f.dndmu=3.0*f.n/2.0/(f.nu-f.m);
272  } else {
273  f.dndmu=3.0*f.n/2.0/f.nu;
274  }
275  f.dsdT=0.0;
276  return;
277  }
278 
279 
280  /** \brief Calculate properties as function of chemical potential
281  at \f$ T=0 \f$
282  */
283  virtual void calc_mu_zerot(fermion_deriv &f) {
284  if (f.non_interacting) { f.nu=f.mu; f.ms=f.m; }
285  if (f.inc_rest_mass) {
286  f.kf=sqrt(2.0*f.ms*(f.nu-f.m));
287  } else {
288  f.kf=sqrt(2.0*f.ms*f.nu);
289  }
290  f.n=f.kf*f.kf*f.kf*f.g/6.0/this->pi2;
291  f.ed=f.g*pow(f.kf,5.0)/20.0/this->pi2/f.ms;
292  if (f.inc_rest_mass) f.ed+=f.n*f.m;
293  f.pr=-f.ed+f.n*f.nu;
294  f.en=0.0;
295 
296  f.dndT=0.0;
297  if (f.inc_rest_mass) {
298  f.dndmu=3.0*f.n/2.0/(f.nu-f.m);
299  } else {
300  f.dndmu=3.0*f.n/2.0/f.nu;
301  }
302  f.dsdT=0.0;
303  return;
304  }
305 
306 
307  /** \brief Calculate properties as function of chemical potential
308  */
309  virtual int calc_mu(fermion_deriv &f, fp_t temper) {
310 
311  if (temper<0.0) {
312  O2SCL_ERR("T<0 in fermion_deriv_nr_tl<fp_t>::calc_mu().",exc_einval);
313  }
314  if (temper==0.0) {
315  calc_mu_zerot(f);
316  return 0;
317  }
318  if (f.non_interacting==true) { f.nu=f.mu; f.ms=f.m; }
319 
320  fp_t pfac2=f.g*pow(f.ms*temper/2.0/this->pi,1.5)/temper, y;
321 
322  if (f.inc_rest_mass) {
323  y=(f.nu-f.m)/temper;
324  } else {
325  y=f.nu/temper;
326  }
327 
328  fp_t half=gsl_sf_fermi_dirac_half(y);
329  fp_t mhalf=gsl_sf_fermi_dirac_mhalf(y);
330  fp_t thalf=gsl_sf_fermi_dirac_3half(y);
331 
332  // Number density:
333  f.n=pfac2*half*temper;
334 
335  // Energy density:
336  f.ed=pfac2*temper*temper*1.5*thalf;
337 
338  if (f.inc_rest_mass) {
339 
340  // Finish energy density
341  f.ed+=f.n*f.m;
342 
343  // entropy density
344  f.en=((f.ed-f.n*f.m)/0.6-(f.nu-f.m)*f.n)/temper;
345 
346  // pressure
347  f.pr=(f.ed-f.n*f.m)/1.5;
348 
349  } else {
350 
351  // entropy density
352  f.en=(5.0*f.ed/3.0-f.nu*f.n)/temper;
353 
354  // pressure
355  f.pr=2.0*f.ed/3.0;
356 
357  }
358 
359  f.dndT=pfac2*(1.5*half-y*mhalf);
360  f.dndmu=pfac2*mhalf;
361  f.dsdT=pfac2*(3.75*thalf-3.0*y*half+y*y*mhalf);
362 
363  return 0;
364  }
365 
366 
367  /** \brief Calculate properties as function of density
368  */
369  virtual int calc_density(fermion_deriv &f, fp_t temper) {
370 
371  if (f.m<0.0 || (f.non_interacting==false && f.ms<0.0)) {
372  O2SCL_ERR2("Mass negative in ",
373  "fermion_deriv_nr_tl<fp_t>::calc_density().",exc_einval);
374  }
375  if (temper<0.0) {
376  O2SCL_ERR2("Temperature less than zero in ",
377  "fermion_deriv_nr_tl<fp_t>::calc_density().",exc_einval);
378  }
379  if (f.n<0.0) {
380  O2SCL_ERR2("Density less than zero in ",
381  "fermion_deriv_nr_tl<fp_t>::calc_density().",exc_einval);
382  }
383 
384  if (temper==0.0) {
386  return 0;
387  }
388  if (f.n==0.0) {
389  if (f.inc_rest_mass) {
390  f.nu=f.m;
391  f.mu=f.m;
392  } else {
393  f.nu=0.0;
394  f.mu=0.0;
395  }
396  f.ed=0.0;
397  f.en=0.0;
398  f.pr=0.0;
399  f.nu=0.0;
400  f.dndT=0.0;
401  f.dndmu=0.0;
402  f.dsdT=0.0;
403  }
404 
405  if (f.non_interacting==true) { f.ms=f.m; f.nu=f.mu; }
406 
407  nu_from_n(f,temper);
408 
409  if (f.non_interacting) { f.mu=f.nu; }
410 
411  calc_mu(f,temper);
412 
413  return 0;
414  }
415 
416 
417  /** \brief Calculate properties with antiparticles as function of
418  chemical potential
419  */
420  virtual int pair_mu(fermion_deriv &f, fp_t temper) {
421 
422  if (f.non_interacting) { f.nu=f.mu; f.ms=f.m; }
423 
424  fermion_deriv antip(f.ms,f.g);
425  f.anti(antip);
426 
427  calc_mu(f,temper);
428 
429  calc_mu(antip,temper);
430 
431  f.n-=antip.n;
432  f.pr+=antip.pr;
433  f.ed+=antip.ed;
434  f.en+=antip.en;
435  f.dsdT+=antip.dsdT;
436  f.dndT+=antip.dndT;
437  f.dndmu+=antip.dndmu;
438 
439  return 0;
440  }
441 
442 
443  /** \brief Calculate properties with antiparticles as function of
444  density
445  */
446  virtual int pair_density(fermion_deriv &f, fp_t temper) {
447  fp_t nex;
448 
449  if (temper<=0.0) {
450  O2SCL_ERR("T=0 not implemented in fermion_deriv_nr().",exc_eunimpl);
451  }
452  if (f.non_interacting==true) { f.nu=f.mu; f.ms=f.m; }
453 
454  nex=f.nu/temper;
455  funct mf=std::bind(std::mem_fn<fp_t(fp_t,fermion_deriv &,fp_t)>
457  this,std::placeholders::_1,std::ref(f),temper);
458 
459  density_root->solve(nex,mf);
460  f.nu=nex*temper;
461 
462  if (f.non_interacting==true) { f.mu=f.nu; }
463 
464  pair_mu(f,temper);
465 
466  return 0;
467  }
468 
469 
470  /// Calculate effective chemical potential from density
471  virtual int nu_from_n(fermion_deriv &f, fp_t temper) {
472 
473  // Use initial value of nu for initial guess
474  fp_t nex;
475  if (f.inc_rest_mass) {
476  nex=-(f.nu-f.m)/temper;
477  } else {
478  nex=f.nu/temper;
479  }
480 
481  // Make a correction if nex is too small and negative
482  // (Note GSL_LOG_DBL_MIN is about -708)
483  if (nex>-GSL_LOG_DBL_MIN*0.9) nex=-GSL_LOG_DBL_MIN/2.0;
484 
485  funct mf=std::bind(std::mem_fn<fp_t(fp_t,fp_t,fp_t)>
487  this,std::placeholders::_1,f.n/f.g,f.ms*temper);
488 
489  // Turn off convergence errors temporarily, since we'll
490  // try again if it fails
491  bool enc=density_root->err_nonconv;
492  density_root->err_nonconv=false;
493  int ret=density_root->solve(nex,mf);
495 
496  if (ret!=0) {
497 
498  // If it failed, try to get a guess from classical_thermo particle
499 
500  classical_thermo cl;
501  cl.calc_density(f,temper);
502  if (f.inc_rest_mass) {
503  nex=-(f.nu-f.m)/temper;
504  } else {
505  nex=-f.nu/temper;
506  }
507  ret=density_root->solve(nex,mf);
508 
509  // If it failed again, add error information
510  if (ret!=0) {
511  O2SCL_ERR("Solver failed in fermion_nonrel::nu_from_n().",ret);
512  }
513  }
514 
515  if (f.inc_rest_mass) {
516  f.nu=-nex*temper+f.m;
517  } else {
518  f.nu=-nex*temper;
519  }
520 
521  return 0;
522  }
523 
524 
525  /** \brief Set the solver for use in calculating the chemical
526  potential from the density */
528  density_root=&rp;
529  return;
530  }
531 
532  /// The default solver for npen_density() and pair_density()
534 
535  /// Return string denoting type ("fermion_deriv_nr")
536  virtual const char *type() { return "fermion_deriv_nr"; };
537 
538  protected:
539 
540 #ifndef DOXYGEN_INTERNAL
541 
542  /// Solver to compute chemical potential from density
544 
545  /// Function to compute chemical potential from density
546  fp_t solve_fun(fp_t x, fp_t nog, fp_t msT) {
547 
548  fp_t nden;
549 
550  // If the argument to gsl_sf_fermi_dirac_half() is less
551  // than GSL_LOG_DBL_MIN (which is about -708), then
552  // an underflow occurs. We just set nden to zero in this
553  // case, as this helps the solver find the right root.
554 
555  if (((-x)<GSL_LOG_DBL_MIN) || !std::isfinite(x)) nden=0.0;
556  else nden=gsl_sf_fermi_dirac_half(-x)*sqrt(this->pi)/2.0;
557 
558  nden*=pow(2.0*msT,1.5)/4.0/this->pi2;
559  return nden/nog-1.0;
560  }
561 
562 
563  /** \brief Function to compute chemical potential from density
564  when antiparticles are included
565  */
566  fp_t pair_fun(fp_t x, fermion_deriv &f, fp_t T) {
567 
568  fp_t nden, y, yy;
569 
570  f.nu=T*x;
571 
572  // 6/6/03 - Should this be included? I think yes!
573  if (f.non_interacting) f.mu=f.nu;
574 
575  if (f.inc_rest_mass) {
576  y=(f.nu-f.m)/T;
577  } else {
578  y=f.nu/T;
579  }
580 
581  nden=gsl_sf_fermi_dirac_half(y)*sqrt(this->pi)/2.0;
582  nden*=f.g*pow(2.0*f.ms*T,1.5)/4.0/this->pi2;
583 
584  yy=nden;
585 
586  if (f.inc_rest_mass) {
587  y=-(f.nu-f.m)/T;
588  } else {
589  y=-f.nu/T;
590  }
591 
592  nden=gsl_sf_fermi_dirac_half(y)*sqrt(this->pi)/2.0;
593  nden*=f.g*pow(2.0*f.ms*T,1.5)/4.0/this->pi2;
594 
595  yy-=nden;
596 
597  yy=yy/f.n-1.0;
598 
599  return yy;
600  }
601 
602 
603 #endif
604 
605  };
606 
607  typedef fermion_deriv_nr_tl<double> fermion_deriv_nr;
608 
609 #ifndef DOXYGEN_NO_O2NS
610 }
611 #endif
612 
613 #endif
o2scl::part_deriv_press_tl::dsdT
fp_t dsdT
Derivative of entropy density with respect to temperature.
Definition: part_deriv.h:135
o2scl::classical_thermo_tl::calc_density
virtual void calc_density(part_tl< fp_t > &p, fp_t temper)
Calculate properties as function of density.
Definition: classical.h:130
o2scl::fermion_deriv_nr_tl::calc_mu
virtual int calc_mu(fermion_deriv &f, fp_t temper)
Calculate properties as function of chemical potential.
Definition: fermion_deriv_nr.h:309
o2scl::fermion_deriv_thermo_tl
Compute properties of a fermion including derivatives [abstract base].
Definition: part_deriv.h:712
o2scl::fermion_deriv_nr_tl
Equation of state for a nonrelativistic fermion.
Definition: fermion_deriv_nr.h:227
o2scl::root::solve
virtual int solve(fp_t &x, func_t &func)=0
o2scl::part_deriv_press_tl::dndT
fp_t dndT
Derivative of number density with respect to temperature.
Definition: part_deriv.h:132
o2scl::part_tl::mu
fp_t mu
Chemical potential.
Definition: part.h:118
o2scl::part_tl::ms
fp_t ms
Effective mass (Dirac unless otherwise specified)
Definition: part.h:122
O2SCL_ERR2
#define O2SCL_ERR2(d, d2, n)
o2scl::fermion_deriv_nr_tl::pair_density
virtual int pair_density(fermion_deriv &f, fp_t temper)
Calculate properties with antiparticles as function of density.
Definition: fermion_deriv_nr.h:446
o2scl::part_tl::n
fp_t n
Number density.
Definition: part.h:112
o2scl::fermion_deriv_nr_tl::calc_density_zerot
virtual void calc_density_zerot(fermion_deriv &f)
Calculate properties as function of density at .
Definition: fermion_deriv_nr.h:255
o2scl::part_tl::non_interacting
bool non_interacting
True if the particle is non-interacting (default true)
Definition: part.h:130
exc_eunimpl
exc_eunimpl
o2scl::part_tl::en
fp_t en
Entropy density.
Definition: part.h:120
o2scl::root
o2scl::fermion_deriv_nr_tl::pair_mu
virtual int pair_mu(fermion_deriv &f, fp_t temper)
Calculate properties with antiparticles as function of chemical potential.
Definition: fermion_deriv_nr.h:420
o2scl::fermion_deriv_thermo_tl< double >::pi
double pi
Desc.
Definition: part_deriv.h:723
o2scl::part_tl::g
fp_t g
Degeneracy (e.g. spin and color if applicable)
Definition: part.h:108
o2scl::part_tl::ed
fp_t ed
Energy density.
Definition: part.h:114
o2scl::fermion_deriv_nr_tl::flimit
fp_t flimit
The limit for the Fermi functions (default 20.0)
Definition: fermion_deriv_nr.h:247
o2scl::fermion_deriv_nr_tl::set_density_root
void set_density_root(root<> &rp)
Set the solver for use in calculating the chemical potential from the density.
Definition: fermion_deriv_nr.h:527
o2scl::fermion_deriv_nr_tl::type
virtual const char * type()
Return string denoting type ("fermion_deriv_nr")
Definition: fermion_deriv_nr.h:536
o2scl::fermion_deriv_nr_tl::nu_from_n
virtual int nu_from_n(fermion_deriv &f, fp_t temper)
Calculate effective chemical potential from density.
Definition: fermion_deriv_nr.h:471
o2scl::fermion_deriv_nr_tl::calc_mu_zerot
virtual void calc_mu_zerot(fermion_deriv &f)
Calculate properties as function of chemical potential at .
Definition: fermion_deriv_nr.h:283
exc_einval
exc_einval
o2scl::fermion_deriv_tl< double >
o2scl::part_tl::inc_rest_mass
bool inc_rest_mass
If true, include the mass in the energy density and chemical potential (default true)
Definition: part.h:128
o2scl::part_deriv_press_tl::dndmu
fp_t dndmu
Derivative of number density with respect to chemical potential.
Definition: part_deriv.h:129
o2scl::fermion_deriv_thermo_tl< double >::pi2
double pi2
Desc.
Definition: part_deriv.h:726
o2scl::fermion_deriv_nr_tl::def_density_root
root_cern def_density_root
The default solver for npen_density() and pair_density()
Definition: fermion_deriv_nr.h:533
o2scl::part_tl::m
fp_t m
Mass.
Definition: part.h:110
o2scl::part_tl::pr
fp_t pr
Pressure.
Definition: part.h:116
o2scl::fermion_tl::kf
fp_t kf
Fermi momentum.
Definition: fermion.h:60
o2scl::part_tl::anti
virtual void anti(part_tl &ax)
Make ap an anti-particle with the same mass and degeneracy.
Definition: part.h:205
o2scl::root_cern
o2scl::fermion_deriv_nr_tl::solve_fun
fp_t solve_fun(fp_t x, fp_t nog, fp_t msT)
Function to compute chemical potential from density.
Definition: fermion_deriv_nr.h:546
o2scl::fermion_deriv_nr_tl::pair_fun
fp_t pair_fun(fp_t x, fermion_deriv &f, fp_t T)
Function to compute chemical potential from density when antiparticles are included.
Definition: fermion_deriv_nr.h:566
o2scl::classical_thermo_tl< double >
O2SCL_ERR
#define O2SCL_ERR(d, n)
o2scl::funct
std::function< double(double)> funct
o2scl::fermion_deriv_nr_tl::calc_density
virtual int calc_density(fermion_deriv &f, fp_t temper)
Calculate properties as function of density.
Definition: fermion_deriv_nr.h:369
o2scl::fermion_deriv_nr_tl::unc
fermion_deriv unc
Storage for the most recently calculated uncertainties.
Definition: fermion_deriv_nr.h:250
o2scl::fermion_deriv_nr_tl::fermion_deriv_nr_tl
fermion_deriv_nr_tl()
Create a fermion with mass m and degeneracy g.
Definition: fermion_deriv_nr.h:232
o2scl::root::err_nonconv
bool err_nonconv
o2scl::part_tl::nu
fp_t nu
Effective chemical potential.
Definition: part.h:124
o2scl::fermion_deriv_nr_tl::density_root
root * density_root
Solver to compute chemical potential from density.
Definition: fermion_deriv_nr.h:536

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).