\documentclass{article}

%other packages
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{physics}

\usepackage[
	style=phys, articletitle=false, biblabel=brackets, chaptertitle=false, pageranges=false, url=true
]{biblatex}

\usepackage{graphicx}
\usepackage{todonotes}
\usepackage{siunitx}

\usepackage[compat=1.0.0]{tikz-feynman}

\usepackage{cleveref}

\title{Notes on the free energy calculation for Owen-Scalapino}
\date{2020-12-20}

\addbibresource{./bibliography.bib}

\graphicspath{{./figures/}}

\newcommand{\pf}{p_{\mathrm{F}}}
\newcommand{\vf}{v_{\mathrm{F}}}
\newcommand{\epsF}{\epsilon_{\mathrm{F}}}
\newcommand{\debye}{\omega_{\mathrm{D}}}
\newcommand{\corr}{\mu^{\ast}}
\newcommand{\dof}{\rho}
\DeclareMathOperator{\sgn}{sgn}

\begin{document}

\maketitle

We are interested in the free energy calculation because we want to be able to correctly limit the range of superconducting behaviour for Owen-Scalapino.
As seen in~\cite{OwenScalapino}, the superconducting state is restricted in its temperature range, with the traditional $T_c$ lowered and with an additional lower critical temperature below which the material returns to the normal state.
There are some subtleties with how we need to handle this calculation in the OS case, which we will discuss.

\section{Superconducting Free energy} \label{sec:scfreeen}

We begin with the standard BCS free energy expression:
\begin{equation}
	F_s = 2 \sum_k \epsilon_k \left(f_k - 2 f_k v_k^2 + v_k^2\right) - \sum_{k, l} V_{kl} u_k v_k u_l v_l (1 - 2 f_k) (1 - 2 f_l) - TS \label{eq:bcs:fs}
\end{equation}

For the distribution function $f_k$, we have
\begin{equation}
	f_k = \frac{1}{1 + \exp\frac{E_k - \corr}{T}},
\end{equation}
where $E_k^2 = \epsilon_k^2 + \Delta^2$.

Our coherence factors are unchanged by the additional $\corr$ constraint in our case, apart from the modified value to $\Delta_k$:
\begin{align}
	\abs{u_k}^2 &= \frac12 \left(1 + \frac{\epsilon_k}{E_k} \right) \\
	\abs{v_k}^2 &= \frac12 \left(1 - \frac{\epsilon_k}{E_k} \right).
\end{align}
These can be obtained most simply from the expression
\begin{equation}
	\tan 2 \theta_k = - \frac{\Delta_k}{\epsilon_k}
\end{equation}
from the December 4 notes, as verification that the coherence factor does not change.
We'll also use
\begin{equation}
	u_k^\ast v_k = \frac12 \frac{\Delta_k}{E_k}
\end{equation}
In \cref{eq:bcs:fs}, incidentally, there are some suppressed magnitude signs, as only one of $u_k$ and $v_k$ can always be made real.
This will be irrelevant here though, I believe, as at best the interaction term in the energy will have the phase cancel out.

We can begin reducing this for code by writing $F_s$ as $A - B - C$ with the three terms in \cref{eq:bcs:fs},
and simplifying each in turn:
\begin{align}
	A &= 2 \sum_k \epsilon_k \left(f_k - 2 f_k v_k^2 + v_k^2 \right) \\
	&= 2 \sum_k \epsilon_k \left(f_k - f_k\left(1 - \frac{\epsilon_k}{E_k} \right) + \frac12 \left(1 - \frac{\epsilon_k}{E_k} \right) \right) \\
	&= 2 \sum_k \epsilon_k \left(f_k\frac{\epsilon_k}{E_k} + \frac12 - \frac12 \frac{\epsilon_k}{E_k} \right) \\
	&= \sum_k \epsilon_k \left(1 + 2 f_k\frac{\epsilon_k}{E_k} - \frac{\epsilon_k}{E_k} \right) \\
\end{align}
Next,
\begin{align}
	B &= \sum_{k, l} V_{kl} u_k v_k u_l v_l (1 - 2 f_k) (1 - 2 f_l) \\
	&= \sum_{k, l} V_{kl} \frac12 \frac{\Delta_k}{E_k} \frac12 \frac{\Delta_l}{E_l} (1 - 2 f_k) (1 - 2 f_l)  \\
	&= \frac14 \sum_{k, l} V_{kl} \frac{\Delta_k}{E_k} \frac{\Delta_l}{E_l} (1 - 2 f_k) (1 - 2 f_l)  \\
\end{align}
Finally, the entropy term:
\begin{align}
	C &= TS \\
	&= -2 T \sum_k f_k \log f_k + (1 - f_k) \log (1 - f_k)
\end{align}

Next, we make our assumptions about the momentum sum.
If we assume that our density of states is constant over the relevant interval, we can replace the integral sums with $\sum_k \rightarrow 2 \int_{0}^{\debye} \dof \dd{\epsilon_k}$.
The $2$ factor here comes from assuming that all of our sums are symmetric in $\epsilon_k$, which is also why our integral begins from $0$.
I don't know if this assumption is correct, but I believe it is (and it handles the intuition that $v_k$ represents our electronlike quasiparticles, which are relevant for $\epsilon_k > \epsF$)
This assumes that we don't consider quasiparticles with energy above $\debye$ as participating in this interaction, and that $\debye \ll \epsF$, which are both reasonable assumptions for this problem.
This means that our term $A$ is then
\begin{align}
	A &= \int_{0}^{\debye} 2 \dof \dd{\epsilon_k} \epsilon_k \left(1 + 2 f_k\frac{\epsilon_k}{E_k} - \frac{\epsilon_k}{E_k} \right) \\
	&= 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} \epsilon_k \left(1 + 2 f_k\frac{\epsilon_k}{E_k} - \frac{\epsilon_k}{E_k} \right) \\
	&= 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} \epsilon_k \left(1 + 2 \left( \frac{1}{1 + \exp\frac{E_k - \corr}{T}} \right)\frac{\epsilon_k}{E_k} - \frac{\epsilon_k}{E_k} \right) \\
	&= \dof\debye^2 + 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} 2 \epsilon_k \left( \frac{1}{1 + \exp\frac{E_k - \corr}{T}} \right)\frac{\epsilon_k}{E_k} - \frac{\epsilon_k^2}{E_k}  \\
	&= \dof\debye^2 + 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} \frac{\epsilon_k^2}{E_k} \left(  \frac{2}{1 + \exp\frac{E_k - \corr}{T}} - 1\right) \\;
	&= \dof\debye^2 + 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} \frac{\epsilon_k^2}{E_k} \frac{2}{1 + \exp\frac{E_k - \corr}{T}} - 2 \dof \int_0^{\debye} \dd{\epsilon_k} \frac{\epsilon_k^2}{\sqrt{k^2 + \Delta^2}} \\
	&= \dof\debye^2 + 2 \dof \int_{0}^{\debye} \dd{\epsilon_k} \frac{\epsilon_k^2}{E_k} \frac{2}{1 + \exp\frac{E_k - \corr}{T}} - \dof \left( \debye \sqrt{\Delta^2 + \debye^2} + \Delta^2 \log \frac{\Delta}{\debye + \sqrt{\Delta^2 + \debye^2}} \right)
\end{align}

Similarly, we can assume that $V_{kl} \rightarrow - V$ and $\Delta_k \rightarrow \Delta$, as part of our standard weak coupling assumption set, giving us
\begin{align}
	B &= V \int_{0}^{\debye} \dof \dd{\epsilon_k} \int_{0}^{\debye} \dof \dd{\epsilon_l} \frac{\Delta_k}{E_k} \frac{\Delta_l}{E_l} (1 - 2 f_k) (1 - 2 f_l)  \\
	&= V \left(\dof \int_{0}^{\debye} \dd{\epsilon_k} \frac{\Delta_k}{E_k} (1 - \frac{2}{1 + \exp\frac{E_k - \corr}{T}}) \right)^2 \\
	&= V \left(\dof \Delta \int_{0}^{\debye} \frac{\dd{\epsilon_k}}{E_k} (1 - \frac{2}{1 + \exp\frac{E_k - \corr}{T}}) \right)^2.
\end{align}
Lastly,
\begin{align}
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \log (1 - f_k) \\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \log (1 - \frac{1}{1 + \exp\frac{E_k - \corr}{T}}) \\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \log (\frac{\exp\frac{E_k - \corr}{T}}{1 + \exp\frac{E_k - \corr}{T}}) \\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \left(\log(\exp\frac{E_k - \corr}{T}) - \log({1 + \exp\frac{E_k - \corr}{T}})\right) \\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \left(\frac{E_k - \corr}{T} - \log({1 + \exp\frac{E_k - \corr}{T}})\right)
\end{align}
We can simplify this slightly further:
\begin{align}
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} f_k \log f_k + (1 - f_k) \left(\frac{E_k - \corr}{T} + \log f_k\right) \\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} \log f_k + (1 - f_k) \frac{E_k - \corr}{T}\\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} -\log({1 + \exp\frac{E_k - \corr}{T}}) + (1 - f_k) \frac{E_k - \corr}{T}\\
	C &= - 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} -\log({\exp\frac{-E_k + \corr}{T}} + 1) - \frac{E_k - \corr}{T} + (1 - f_k) \frac{E_k - \corr}{T}\\
	C &= 4 T \int_{0}^{\debye} \dof \dd{\epsilon_k} \log({\exp\frac{-E_k + \corr}{T}} + 1) + f_k \frac{E_k - \corr}{T}\\
\end{align}

\printbibliography

\end{document}
