Integration of Interpolated function take an unacceptable amount of time
The code in the question takes about 40 seconds on my computer. Turning off SymbolicProcessing
, as suggested by Tim Laska in a comment, reduces the run time to about 16 seconds. Although reducing MaxRecursion
further reduces the run time, it also generates error messages that the integration is not converging to sufficient accuracy. (Whether this matters from a practical point of view depends on how the results will be used.)
Instead, let us take advantage of the fact that NIntegrate
in Plot
repeatedly performs the same integral over the range {2.66*^17, 3.78*^18}
to break the integral into two pieces.
int = NIntegrate[ElectronNumberDensity[eta], {eta, 2.66*^17, 3.78*^18}]
Plot[int + NIntegrate[ElectronNumberDensity[eta], {eta, η, 2.66*^17},
Method -> {Automatic, "SymbolicProcessing" -> 0}], {η, 1.47*^17, 2.66*^17},
ImageSize -> Large, LabelStyle -> {Bold, Black, 15}]
which produces the desired plot in under 5 seconds.
The plot agrees well with that produced by the code in the question.
Addendum
Even faster is to use NDSolveValue
instead of NIntegrate
.
NDSolveValue[{g'[eta] == -ElectronNumberDensity[eta], g[2.66*^17] == int},
g, {eta, 1.47*^17, 2.66*^17}]
Plot[%[η], {η, 1.47*^17, 2.66*^17}, ImageSize -> Large, LabelStyle -> {Bold, Black, 15}]
which produces the same curve in less than a second.
This
c = 2.99792*10^5;
A = 3.87624*10^-14;
FreeElectronFractionData = {...YourData...};
Clear[\[Eta]]; (* just in case there was a prior assignment to Eta *)
redShift = 6.64*^18^2/((c - Sqrt[c]*Sqrt[c - 2.*A*\[Eta]])/A)^2 - 1.;
FreeElectronFraction = Interpolation[FreeElectronFractionData, InterpolationOrder -> 1];
ElectronNumberDensity[\[Eta]_] := FreeElectronFraction[redShift ]*1.42*^-7*(1. +redShift^3);
Plot[NIntegrate[ElectronNumberDensity[eta], {eta, \[Eta], 3.78*^18}, MaxRecursion -> 15], {\[Eta], 1.47*^17, 2.66*^17}]
ListPlot[FreeElectronFractionData]
completes in a less than ten seconds, mostly I believe because I changed
FreeElectronFraction := Interpolation[...]
to
FreeElectronFraction = Interpolation[...]
which means the calculation of the interpolating function will be done only once and the result will then be repeatedly used instead of calculating the interpolation function thousands of times.
I also moved the location of your calculation of redShift
which may or may not have been necessary, but I was getting recursion errors and fixed those by moving the redShift
calculation.
Please check all this very carefully to make certain the result is correct.