Compiling the VoigtDistribution PDF

Here is a compiled implementation of the Voigt profile function, based on an approximation derived by Chiarella and Reichel and improved by Abrarov, Quine and Jagpal:

voigt = With[{n = 24, τ = 12},
             With[{d = N[Range[n] π/τ], b = N[Exp[-(Range[n] π/τ)^2]],
                   s = N[PadRight[{}, n, {-1, 1}]],
                   sq = N[Sqrt[2]], sp = N[Sqrt[2 π]]},
                  Compile[{{δ, _Real}, {σ, _Real}, {x, _Real}},
                          Module[{z = (x + I δ)/(σ sq), e}, e = Exp[I τ z];
                                 Re[(I (1 - e)/(τ z) + (2 I z/τ)
                                    b.((e s - 1)/((d + z) (d - z))))]/(σ sp)],
                          RuntimeAttributes -> {Listable}]]];

The compiled function performs remarkably well on my box:

vp = PDF[VoigtDistribution[1, 3/2]];
AbsoluteTiming[v1 = vp[N[Range[-7, 7, 1/50]]];]
   {3.70813, Null}

AbsoluteTiming[v2 = voigt[1, 3/2, Range[-7, 7, 1/50]];]
   {0.053813, Null}

Norm[v1 - v2, ∞]
   1.38778*10^-16

Further speed could probably be squeezed out by explicitly expanding out the complex arithmetic into its real and imaginary components, but the resulting code will look comparatively unwieldy.


From the documentation, we have a pseudo-Voigt Distribution that might be used as an approximation. This might be useful as a basis for making a compilable function.

PseudoVoigtDistribution[δ_, σ_] := 
 Block[{g = (δ^5 + σ^5 + 2.69296 σ^4 δ + 2.42843 σ^3 δ^2 + 
        4.47163 σ^2 δ^3 + 0.07842 σ δ^4)^(1/5), η},
       η = δ/g;
       η = η*(1.36603 - 0.47719 η + 0.11116 η^2);
       MixtureDistribution[{1 - η, η}, {NormalDistribution[0, g], 
                           CauchyDistribution[0, g]}]
       ]