Computation of a general integral of Fermi‐Dirac distribution by McDougall‐Stoner method

We extended the method of McDougall and Stoner (1938) [15] to the computation of a general integral of the Fermi–Dirac distribution, $F(\eta )$ . When $\eta >0$ , the new method splits $F(\eta )$ into a sum of three parts, $A(\eta )\text{,}B(\eta )$ , and $C(\eta )$ , and integrates them exactly and/or numerically by means of the double‐exponential quadrature rules. As η increases, $B(\eta )/F(\eta )$ damps algebraically while $C(\eta )/F(\eta )$ decays exponentially. Thus, they can be ignored when η exceeds certain threshold values depending on the input error tolerance. When $A(\eta )$ is exactly computable and η is sufficiently large such that $B(\eta )$ and $C(\eta )$ are negligible, the new method was 60–1000 times faster than the numerical integration of $F(\eta )$ . Even if $B(\eta )$ and/or $C(\eta )$ are not negligible, their smallness and rapid decay significantly accelerate their numerical quadrature so that the speed‐up factor becomes 3–20 unless η is less than 2.0–4.5. On the other hand, when $A(\eta )$ is numerically integrated, the speed‐up factor diminishes to 2–5. Still, the superiority of the new method to the numerical integration of $F(\eta )$ is unchanged.