# Numerically convolve lists in mathematica

I have a function $$\chi(x) = \frac{a^2}{a^2 + \left(\frac{d^{2}}{x} + x\left( \frac{ab}{x^{2} - c^2} - 1 \right) \right)^{2}}$$ which I would like to convolve with the Fourier transform of the Hanning window, $$H(t)\rightarrow\tilde{H}(\omega)$$. I have attempted to do this 'analytically' in Mathematica with

Convolve[\chi[x], FourierTransform[HanningWindow[t], t, x], x, y]


However this never completes even with assumptions, such as Assumptions-> x > 0 && y > 0 && a > 0 && b > 0 && c > 0. Maybe this expression is too difficult, there is maybe no closed form or maybe I don't know enough about how to use Convolve[...] effectively.

If anyone has any thoughts about this problem that would be very helpful. Here is the function in, $$\chi(x)$$, Mathematica code:

\[Chi][a_, b_, c_, d_, x_] := a^2/(
a^2 + (d^2/x + x ((a b)/(x^2 - c^2) - 1))^2)


My next approach is to do this numerically. So I will generate two lists one for $$\chi(x)$$, where I define values for $$a$$, $$b$$, and $$c$$, and one list for $$\tilde{H}(\omega)$$. I choose some arbitary resolution in $$\omega$$ as $$\delta \omega$$.

I then want to make the convolution of the two lists $$\{\chi(\omega_{i})\} * \{\tilde{H}(\omega_{i})\}$$, and maybe this will give me the information I am interested in.

I haven't been able to discover how to make such a convolution with lists in Mathematica -- is this possible?

• Check ListConvolve. – Daniel Lichtblau Aug 14 '19 at 15:31
• Yours code is wrong? Probably you want:FourierTransform[HannWindow[x], x, ω]? – Mariusz Iwaniuk Aug 14 '19 at 16:06
• @MariuszIwaniuk yes that is true, that's just a typo though. The convolution still doesn't complete. – Q.P. Aug 14 '19 at 16:07
• I absolutely did! Thanks, I have corrected. – Q.P. Aug 14 '19 at 16:17
• Can you show your function definition for chi? – bill s Aug 14 '19 at 19:16

Both your HannWindowand $$\chi$$ function are Indeterminate at $$x=0$$ so I created a function that uses the Limit at $$x=0$$ or I rearranged the function to eliminate the divide by zero.

I don't know if it is necessary, but I seem to have had better success if I make the function cyclical by mirroring it about it's endpoint before doing the convolution.

The following code will create a numerical ListConvolve and a Manipulate so you can see the effect of parameter changes.

window[ω_][x_] =
If[x == 0, 1/(
2 Sqrt[2 π] ω), (-((
2 I E^(-((I x)/(
2 ω))) (-1 + E^((
I x)/ω)) π^2 ω^2 (UnitStep[-ω] -
UnitStep[ω]))/(x^3 - 4 π^2 x ω^2)) +
2 π DiracDelta[x] UnitStep[-ω] UnitStep[ω])/
Sqrt[2 π]];
χ[a_, b_, c_, d_][
x_] := (a^2*x^2)/(d^4 + x^4 + (2*a*b*x^4)/(c^2 - x^2) +
a^2*(x^2 + (b^2*x^4)/(c^2 - x^2)^2) +
2*d^2*x^2*(-1 + (a*b)/(-c^2 + x^2)));
delta[x0_, xf_, n_] := (xf - x0)/(n - 1)
grid[x0_, xf_, n_] := N@Range[x0, xf, delta[x0, xf, n]]
cyclic = N@(#~Join~Reverse[#] &)@(dummy /@ #) &;
discreteConvolve[startx_, finalx_, npoints_, w_, m_] :=
Module[{gr, g, wdigitized, mdigitized, conv, discreteConv},
(* Create the cyclic grid *)
gr = grid[startx, finalx, npoints];
g = cyclic@gr;
(* Digitize the window/kernel *)
wdigitized = w @@@ g;
(* Normalization *)
wdigitized = wdigitized/Total[wdigitized];
(* Digitize the model *)
mdigitized = m @@@ g;
(* Perform the convolution *)
conv = ListConvolve[wdigitized, mdigitized, {1, -1}, 0];
discreteConv = Transpose[{gr, conv~Take~npoints}];
(* Display the results *)
discreteConv]
Manipulate[
Show[ListPlot[
discreteConvolve[0, finalx, npoints,
window[ω], χ[a, b, c, d]],
PlotRange -> {{0, finalx}, {0, 1}},
PlotLegends -> {"Convolution"}],
Plot[χ[a, b, c, d][x], {x, 0, finalx},
PlotRange -> {{0, finalx}, {0, 1}}, PlotStyle -> {Red, Dashed},
PlotLegends -> {"Model"}]],
{{ω, 1}, 0.01, 1},
{{a, 1}, 0, 6},
{{b, 1}, 0, 6},
{{c, 1}, 0, 6},
{{d, 1}, 0, 6},
{{finalx, 30}, 10, 100},
{{npoints, 200}, 100, 5000, 100}
]


• Wow this is really neat!! It's okay though as $x$ is always real and positive. Seriously this is really nice. – Q.P. Aug 15 '19 at 17:01
• If I could ask one question before I accept your answer, what is the meaning of the parameter omega you introduce. I.e. I guess this is related to the Fourier transform if the window, but why is it a variable? Is this the size of the transformed window in frequency space? – Q.P. Aug 16 '19 at 12:21
• @QuantumPenguin Thank you for your kind words. You are correct. Omega was a term that I introduced to vary the width of the window function. I was not sure what the nominal values of the a-d parameters were and with my range choice, an omega of 1 excessively broadened the peaks. I did it to convince myself that ListConvolve would sharpen as I narrowed the window function. – Tim Laska Aug 16 '19 at 12:38