2
$\begingroup$

I'd like to illustrate the mathematical process of convolution ($f(x) \star g(x)$) with an interactive Manipulate, but the evaluation is far too slow for a class lecture.

Here's the basic (unoptimized) code:

Manipulate[
 GraphicsColumn[{
   Plot[Evaluate@Convolve[f[x], g[x], x, y], {y, 0, 10},
    AxesLabel -> {"z", " (f\[Star]g)(z)"},
    PlotRange -> {-2, 2}, Exclusions -> None,
    Epilog -> {If[(hh = 
          Evaluate@Convolve[f[x], g[x], x, y] /. y -> z) > 0, Green, 
       Red], PointSize[0.04], Point[{z, hh}]}],
   Plot[f[x] g[z - x], {x, 0, 10}, 
    AxesLabel -> {"x", " f(x) g(z-x)"}, PlotRange -> {-1, 1}, 
    Exclusions -> None, Filling -> Axis, 
    FillingStyle -> {Directive[Opacity[0.3], Red], 
      Directive[Opacity[0.3], Green]}],
   Plot[g[z - x], {x, 0, 10}, AxesLabel -> {"x", "g(z-x)"}, 
    PlotRange -> {-1, 1}, 
    Epilog -> {Black, Line[{{z, 0}, {z, -.2}}], 
      Text[Style["z", Italic, 24, FontFamily -> "Times"], {z, -.33}]},
     Exclusions -> None],
   Plot[f[x], {x, 0, 10}, PlotRange -> {-1, 1}, 
    AxesLabel -> {"x", " f(x)"}, Exclusions -> None]}],
 {f, {Sin, SquareWave, TriangleWave}},
 {g, {UnitBox, HeavisideLambda[#/3] &, UnitTriangle[#/5] &}},
 {{z, 0}, 0, 10},
 ContinuousAction -> False]

It is simplest to read the four-part figure from the bottom. The bottom plot is $f(x)$. The next-to-bottom is $g(z-x)$, where $z$ is controllable by a slider. The second plot from the top is the product $f(x) g(z-x)$. I've used shading to highlight the positive values (green) and negative values (red). The convolution is the full integral of this function, i.e., the green shaded area minus the red shaded area.

For each chosen value of $z$, the value of that integral is plotted as a point on the top plot, again, colored green if positive, and red if negative. The curve plotted in the top plot is the convolution, $(f(x) \star g(x))[z]$.

This code is acceptably fast when initialized to $f = \sin x$ and $g = $ UnitBox, as you can confirm by running the code. I've tried to speed things by using ContinuousAction -> False, which helps just a bit.

But the execution is exceedingly slow for more general functions.

I believe that the convolution is re-computed for each setting of $z$, which of course is highly inefficient.

There must be some way in which the convolutions are pre-computed and stored, or numerically evaluated, or otherwise processed to speed the overall evaluation, but I haven't found a way to make the code interactive in real time.

$\endgroup$
2
  • 2
    $\begingroup$ I have not tried all combinations, but a number of your $f$ and $g$ do not convolve. For example Convolve[Sin[x],HeavisideLambda[x/3],x,y] and Convolve[TriangleWave[x], UnitTriangle[x/5], x, y]. Mathematica simply returns them unevaluated. So no point of passing the result to Plot command, as it will hang. It will be better to first test these functions outside of Manipulate to make sure they work OK and fast, then add them to Manipulate. If slowness then result, then we know it is because of Manipulate setting. $\endgroup$
    – Nasser
    Commented Dec 10, 2023 at 1:37
  • $\begingroup$ There is an interactive example on the doc page. This page has a few responsive interactive examples as well. $\endgroup$
    – Syed
    Commented Dec 10, 2023 at 2:50

1 Answer 1

7
$\begingroup$

As pointed by Nasser in comment above, Convolve is having difficulty in handling certain combination of f and g, for example:

Convolve[Sin[x], HeavisideLambda[x/3], x, y]
(* Input returned *)

Though not documented, when y is a numeric value, Convolve can compute the convolution to some degree:

Convolve[Sin[x], HeavisideLambda[x/3], x, 1.2] // AbsoluteTiming
(* {0.692999, 1.2365} *)

Convolve[Sin[x], HeavisideLambda[x/3], x, 3] // AbsoluteTiming
(* {0.590948, 4/3 Sin[3/2]^2 Sin[3]} *)

But this is just too slow. What's more, this hidden functionality is unstable:

Convolve[Sin[x], HeavisideLambda[x/3], x, 3.]
(* Input returned *)

So we need a better numeric solver. It's straightforward to define numeric convolution via NIntegrate:

Clear[NConvolve]; 
NConvolve[f_, g_, x_, y_?NumericQ] := 
 NConvolve[f, g, x, y] = 
  NIntegrate[f (g /. x -> y - x), {x, -∞, ∞}, AccuracyGoal -> 5]

The AccuracyGoal -> 5 is added to avoid the difficulty caused by zero integration. To speed things up, we pre-compute those convolutions:

{L, R} = {0, 10}; step = 1/10;

rule = HeavisideLambda -> UnitTriangle;

Do[
  NConvolve[f[x], g[x] /. rule, x, 
   y], {f, {Sin, SquareWave, TriangleWave}}, {g, {UnitBox, HeavisideLambda[#/3] &, 
    UnitTriangle[#/5] &}}, {y, L, R, step}] // AbsoluteTiming
(* {53.0485, Null} *)

The rule is added because HeavisideLambda is not suitable for numeric calculation.

Remark

The convolution only needs to be calculated once. We can then store the result using e.g.

DumpSave["convolve.mx", NConvolve]

and load the result using e.g.

<<convolve.mx

The next optimization is to use Grid instead of GraphicsColumn, which is slow:

Clear[convplot, fgplot, gplot, fplot]
convplot[f_, g_] :=(*convplot[f,g]=*)
 ListLinePlot[Table[{y, NConvolve[f[x], g[x] /. rule, x, y]}, {y, L, R, step}], 
  AxesLabel -> {"z", " (f⋆g)(z)"}, PlotRange -> {-2, 2}]
point[f_, g_, z_] := 
 Graphics@{If[(hh = NConvolve[f[x], g[x] /. rule, x, z]) > 0, Green, Red], 
   PointSize[0.04], Point[{z, hh}]}
fgplot[f_, g_, z_] :=(*fgplot[f,g,z]=*)
 Plot[f[x]  g[z - x], {x, L, R}, AxesLabel -> {"x", " f(x) g(z-x)"}, 
  PlotRange -> {-1, 1}, Exclusions -> None, Filling -> Axis, 
  FillingStyle -> {Directive[Opacity[0.3], Red], Directive[Opacity[0.3], Green]}]
gplot[g_, z_] :=(*gplot[g,z]=*)
 Plot[g[z - x], {x, L, R}, AxesLabel -> {"x", "g(z-x)"}, PlotRange -> {-1, 1}, 
  Epilog -> {Black, Line[{{z, 0}, {z, -.2}}], 
    Text[Style["z", Italic, 24, FontFamily -> "Times"], {z, -.33}]}, 
  Exclusions -> None]
fplot[f_] :=(*fplot[f]=*)
 Plot[f[x], {x, L, R}, PlotRange -> {-1, 1}, AxesLabel -> {"x", " f(x)"}, 
  Exclusions -> None]

Manipulate[
 Grid[{{Show[convplot[f, g], point[f, g, z]], fgplot[f, g, z]}, {gplot[g, z], 
    fplot[f]}}], {f, {Sin, SquareWave, TriangleWave}}, {g, {UnitBox, 
   HeavisideLambda[#/3] &, UnitTriangle[#/5] &}}, {z, L, R, step}]

Notice I've set an explicit step in Manipulate so the pre-computed values will be called. We can add memorization for convplot, etc to speed up the code further, but it's not quite necessary in my view, because the current implementation is already smooth even on my old laptop:

enter image description here

$\endgroup$
1
  • $\begingroup$ Oh how very nice. Thanks for all your efforts. A hearty $+1$ and $\checkmark$. I'd been working on pre computing convolutions and then merely reading off/implementing the function. (My students will love this.) $\endgroup$ Commented Dec 10, 2023 at 18:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.