Adds free energy case to mathematica package

src/wl/owenScalapinoFreeEnergy.wl

@@ -6,6 +6,8 @@ entropy::usage = "entropy[delta, temp, mustar, n0] returns superconducting entro
 ATerm::usage = "Calculates A term from notes";
 BTerm::usage = "B term";
 CTerm::usage = "C term";
+DTerm::usage = "D term";
+Fs::usage "Free energy for sc case";
 
 Begin["`Private`"];
 
@@ -51,7 +53,19 @@ BTerm[delta_?NumericQ, temp_?NumericQ, mustar_?NumericQ,
 			BTermIntegrand[en, delta, temp, mustar], {en, 0, debyeFreq}] )^2;
 
 CTerm[delta_?NumericQ, temp_?NumericQ, mustar_?NumericQ,
-  N0_?NumericQ] := -entropy[delta, temp, mustar, N0] * temp;
+	N0_?NumericQ] := -entropy[delta, temp, mustar, N0] * temp;
+
+DTerm[delta_?NumericQ, temp_?NumericQ, mustar_?NumericQ,
+	debyeFreq_?NumericQ, N0_?NumericQ] := - 4 * mustar * N0 *
+		NIntegrate[fs[en, delta, temp, mustar], {en, 0, debyeFreq}];
+
+(* Bring it home *)
+Fs[delta_?NumericQ, temp_?NumericQ, mustar_?NumericQ,
+	debyeFreq_?NumericQ, N0_?NumericQ, V0_?NumericQ] :=
+		ATerm[delta, temp, mustar, debyeFreq] +
+				BTerm[delta, temp, mustar, debyeFreq, N0, V0] +
+				CTerm[delta, temp, mustar, N0] +
+				DTerm[delta, temp, mustar, debyeFreq, N0];
 
 End[]; (* `Private` *)
 

src/wltest/BasicTest.wlt

@@ -35,4 +35,18 @@ VerificationTest[
 	, TestID -> "C term"
 ];
 
+VerificationTest[
+	DTerm[1.3, .5, .2, 100, 1]
+	,
+	-0.09439746161647347
+	, TestID -> "D term"
+];
+
+VerificationTest[
+	Fs[.5, .5, .247777, 50, 1, .2]
+	,
+	1.7128117685372855
+	, TestID -> "Fs Test"
+];
+
 EndTestSection[];