From 65bdf79a8fd5d039b28a6490ed365288c5bb8e00 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Tue, 1 Sep 2026 23:51:54 -0400 Subject: [PATCH 1/5] added notebook to explain math ffor convolution --- docs/overlap_add_math.ipynb | 222 ++++++++++++++++++++++++++++++++++++ 1 file changed, 222 insertions(+) create mode 100644 docs/overlap_add_math.ipynb diff --git a/docs/overlap_add_math.ipynb b/docs/overlap_add_math.ipynb new file mode 100644 index 0000000..346a29e --- /dev/null +++ b/docs/overlap_add_math.ipynb @@ -0,0 +1,222 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "e984e902-bf76-43ae-a5d5-55e59357b6bf", + "metadata": {}, + "source": [ + "# From Linear Convolution to Overlap-Add Method!" + ] + }, + { + "cell_type": "markdown", + "id": "e9725572-a39b-4fb5-992c-5f3d22a80797", + "metadata": {}, + "source": [ + "This notebook explains the math behind the overlap-add method in convolution. The content is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/)" + ] + }, + { + "cell_type": "markdown", + "id": "03f9d2f8-1560-4a5d-93a2-ec82fb422f6f", + "metadata": {}, + "source": [ + "Most textbooks use $x$ and $h$ notation in the convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. We would also like to introduce $Q'$, which is simply the reverse of $Q$." + ] + }, + { + "cell_type": "markdown", + "id": "3dda604a-3a3e-42f4-90c7-7b01ada16337", + "metadata": {}, + "source": [ + "# Linear Convolution of Two Finite-Length Sequences" + ] + }, + { + "cell_type": "markdown", + "id": "9ee474bd-1229-4749-80fe-f90c4a135931", + "metadata": {}, + "source": [ + "Consider two finte-length sequences $Q'$ (of length $m$) and $T$ (of length $n$). Their linear convolution, $C$, can be computed as follows:" + ] + }, + { + "cell_type": "markdown", + "id": "ad4a63b5-2d91-4875-b4c2-bc7ce11b6af9", + "metadata": {}, + "source": [ + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$" + ] + }, + { + "cell_type": "markdown", + "id": "8028ee30-3244-4c79-86ca-22e2b8876694", + "metadata": {}, + "source": [ + "where $T[.]$ and $Q'[.]$ are zeros for any index that is outside of their range. So, $T[i]$ is zero if $i$ is outside of the range $0 \\le i \\le n-1$. Hence, the equation above is equivalent to:\n", + "\n", + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]} $$" + ] + }, + { + "cell_type": "markdown", + "id": "593a8276-49cb-43d2-9300-8d29e67d3104", + "metadata": {}, + "source": [ + "Note that $Q'[idx-i]$ is not zero when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. \n", + "\n", + "Therefore:\n", + "\n", + "* $idx \\ge i, \\quad i \\ge 0 \\implies idx \\ge 0$\n", + "* $idx \\le i+m-1, \\quad i \\le n-1 \\implies idx \\le n+m-2$" + ] + }, + { + "cell_type": "markdown", + "id": "e638d3d7-cf2a-4355-aff3-a9192400defb", + "metadata": {}, + "source": [ + "This shows that the linear convolution $C$ has values at indices $\\set{0, 1, ..., n + m - 2}$, and it is 0 otherwise. Therefore, the length of output in linear convolution is $n+m-1$. \n", + "\n", + "Furthermore, let's check out the values for different indices:\n", + "* $idx=0 \\implies C[0]=T[0]Q'[0]$\n", + "* $idx=1 \\implies C[0]=T[0]Q'[1] + T[1]Q'[0]$\n", + "* ...\n", + "* $idx=m-1 \\implies C[m-1]=T[0]Q'[m-1] + T[1]Q'[m-2] + ... + T[m-1]Q'[0] = T_{0}.Q$\n", + "* $idx=m \\implies C[m]=T[1]Q'[m-1] + T[2]Q'[m-2] + ... + T[m]Q'[0] = T_{1}.Q$\n", + "* ...\n", + "* $idx=n-1 \\implies C[n-1]=T[n-m]Q'[(n-1)-(n-m)] + T[n-m+1]Q'[(n-1)-(n-m+1)] + ... + T[n-1]Q'[0] = T_{n-m+1}.Q$\n", + "* ...\n", + "* $idx=n+m-2 \\implies C[n+m-2]=T[n-1]Q'[(n+m-2)-(n-1)] = T[n-1]Q'[m-1]$\n", + "\n", + "As observed, when $m-1 \\le idx \\le n-1$, $C[idx]$ becomes the dot product between a subsequnce of $T$ and $Q$, which is the reverse of $Q'$. Therefore, the `range(m-1,n)` gives the sliding dot product between $Q$ and $T$." + ] + }, + { + "cell_type": "markdown", + "id": "48c7f9d8-1b7c-4abc-886b-27a015b78ed6", + "metadata": {}, + "source": [ + "**How can we leverage this to speed up the computation of sliding dot product?**" + ] + }, + { + "cell_type": "markdown", + "id": "52bcf566-3ac0-483d-81a3-57fa3b5766c3", + "metadata": {}, + "source": [ + "# Option I: Convert to Circular Convolution and use FFT-IFFT" + ] + }, + { + "cell_type": "markdown", + "id": "2c4176e2-60f3-410f-ad2d-6c1d29cca6cf", + "metadata": {}, + "source": [ + "Circular convolution is defined between two sequences that are both periodic and their period are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n", + "\n", + "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{Q'}[i] \\times \\tilde{T}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", + "\n", + "where, $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$. $C_{N}$ represents one period of N-Circular convolution. This can also be computed via FFT, i.e. $$C_{N} = IFFT( FFT(\\tilde{Q'}_{N}) \\times FFT(\\tilde{T}_{N}) ) $$" + ] + }, + { + "cell_type": "markdown", + "id": "e10c5304-a1a1-4cbd-9b26-266e229bdc57", + "metadata": {}, + "source": [ + "If there is a way to compute the linear convolution via circular convolution, then we can take advantage of FFT by using eq (4). The good news is that there is a way! The linear convolution between $Q'$ and $T$ can be obtained by performing N-circular convolution between $\\tilde{Q'}$ and $\\tilde{T}$, where:\n", + "\n", + "* $N \\ge n + m - 2$\n", + "* $\\tilde{Q'}_{N}$ is $Q'$ but zero-padded with $N-m$ zeros\n", + "* $\\tilde{T}_{N}$ is $T$ but zero-padded with $N-n$ zeros" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "b6529137-fcdc-4d65-aa93-da3746739f94", + "metadata": {}, + "outputs": [], + "source": [ + "# WIP" + ] + }, + { + "cell_type": "markdown", + "id": "97fcc328-35f9-4b66-9598-1b5d8f1f97ac", + "metadata": {}, + "source": [ + "# Option II: Overlap-add method\n", + "This is a divide and conquer algorithm. It applies divide on \"linear convolution\" and conquer each via circular convolution (See: Option I)" + ] + }, + { + "cell_type": "markdown", + "id": "1e45ccf7-1564-4c41-b32e-4de9b1b201a3", + "metadata": {}, + "source": [ + "## Linearity in Linear Convolution" + ] + }, + { + "cell_type": "markdown", + "id": "a084578e-cb20-437c-a70c-3e33f8da8ac3", + "metadata": {}, + "source": [ + "Suppose the array $T$ can be written as $T_{1} + T_{2}$. In other words: $T[i] = T_{1}[i] + T_{2}[i]$, then:\n", + "\n", + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$\n", + "\n", + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{(T_{1}[i] \\times Q'[idx-i] + T_{2}[i] \\times Q'[idx-i])} $$\n", + "\n", + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T_{1}[i] \\times Q'[idx-i]} + \\sum_{i=-\\infty}^{\\infty}{T_{2}[i] \\times Q'[idx-i]}$$\n", + "\n", + "$$ C[idx] = C_{1}[idx] + C_{2}[idx] $$" + ] + }, + { + "cell_type": "markdown", + "id": "1f5be13b-82a7-46d9-9cb1-229c878e67f8", + "metadata": {}, + "source": [ + "This shows that a linear convolution has the property of \"linearity\". Now, we use this property to show that we can compute the linear convolution between $Q'$ and long $T$ by breaking $T$ into smaller parts." + ] + }, + { + "cell_type": "markdown", + "id": "fecee5fc-3c10-4c2e-8c66-bce17c3650bd", + "metadata": {}, + "source": [] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "acf5516d-e52d-457c-a8c5-f4d1b7ae8f1c", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.14.4" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} From f10007723333fb1133840f0aa122a315360a635e Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Thu, 3 Sep 2026 01:30:28 -0400 Subject: [PATCH 2/5] complete Option 1 part and added bonus --- docs/overlap_add_math.ipynb | 96 +++++++++++++++++++++++++++++++++---- 1 file changed, 88 insertions(+), 8 deletions(-) diff --git a/docs/overlap_add_math.ipynb b/docs/overlap_add_math.ipynb index 346a29e..8d0b9e7 100644 --- a/docs/overlap_add_math.ipynb +++ b/docs/overlap_add_math.ipynb @@ -115,9 +115,9 @@ "source": [ "Circular convolution is defined between two sequences that are both periodic and their period are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n", "\n", - "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{Q'}[i] \\times \\tilde{T}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", + "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", "\n", - "where, $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$. $C_{N}$ represents one period of N-Circular convolution. This can also be computed via FFT, i.e. $$C_{N} = IFFT( FFT(\\tilde{Q'}_{N}) \\times FFT(\\tilde{T}_{N}) ) $$" + "where, $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$. $C_{N}$ represents one period of N-Circular convolution. This can also be computed via FFT, i.e. $C_{N} = IFFT( FFT(\\tilde{T}_{N}) \\times FFT(\\tilde{Q'}_{N}) )$." ] }, { @@ -125,21 +125,101 @@ "id": "e10c5304-a1a1-4cbd-9b26-266e229bdc57", "metadata": {}, "source": [ - "If there is a way to compute the linear convolution via circular convolution, then we can take advantage of FFT by using eq (4). The good news is that there is a way! The linear convolution between $Q'$ and $T$ can be obtained by performing N-circular convolution between $\\tilde{Q'}$ and $\\tilde{T}$, where:\n", + "If there is a way to compute the linear convolution via circular convolution, then we can use FFT to compute the linear convolution. Actually, there is a way! The linear convolution between $Q'$ and $T$ is equivalent to the N-circular convolution between $\\tilde{Q'}$ and $\\tilde{T}$, where:\n", "\n", - "* $N \\ge n + m - 2$\n", + "* $N \\ge n + m - 1$, to cover the $n+m-1$ data points of linear convolution\n", "* $\\tilde{Q'}_{N}$ is $Q'$ but zero-padded with $N-m$ zeros\n", "* $\\tilde{T}_{N}$ is $T$ but zero-padded with $N-n$ zeros" ] }, { - "cell_type": "code", - "execution_count": 2, + "cell_type": "markdown", "id": "b6529137-fcdc-4d65-aa93-da3746739f94", "metadata": {}, - "outputs": [], "source": [ - "# WIP" + "So:\n", + "\n", + "$$\n", + "\\tilde{Q'}_{N}[i] = \\begin{cases} \n", + " Q'[i] & \\text{if } 0 \\le i \\le m-1 \\\\\n", + " 0 & \\text{if } m \\le i \\le N-1\n", + " \\end{cases}\n", + "$$\n", + "\n", + "$$\n", + "\\tilde{T}_{N}[i] = \\begin{cases} \n", + " T[i] & \\text{if } 0 \\le i \\le n-1 \\\\\n", + " 0 & \\text{if } n \\le i \\le N-1\n", + " \\end{cases}\n", + "$$" + ] + }, + { + "cell_type": "markdown", + "id": "5ce64782-66d0-413b-a116-9ac3555966b6", + "metadata": {}, + "source": [ + "Let's break down the circular convolution equation:\n", + "\n", + "$$C_{N}[idx] = \\sum_{i=0}^{m-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=m}^{n-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=n}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", + "\n", + "$$C_{N}[idx] = \\sum_{i=0}^{m-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=m}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", + "\n", + "$$C_{N}[idx] = \\sum_{i=0}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$" + ] + }, + { + "cell_type": "markdown", + "id": "a7e66371-4cd4-4673-957c-b7c548dcbfde", + "metadata": {}, + "source": [ + "Let's compare the linear convolution in eq(2) and the N-circular convolution in eq(9). To show that they give the same value, we need to compare $Q'[idx-i]$ with $\\tilde{Q'}[(idx-i)_{N}]$.\n", + "\n", + "Let $j$ denote $idx-i$. So, all we need to do is to compare $Q'[j]$ with $\\tilde{Q'}[(j)_{N}]$.\n", + "\n", + "**Case A: $j \\ge 0$**
\n", + "\n", + "In this case, $max(j) = max(idx) - min(i) = (N-1) - 0 = N-1 \\implies \\tilde{Q'}[(j)_{N}]=Q'_{N}[j]$\n", + "\n", + "\n", + "**Case B: $j < 0$**
\n", + "\n", + "In this case, $Q'[j]=0$.\n", + "\n", + "In this case, $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$. Note that: \n", + "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge (0 - (n-1)) + N$$ \n", + "\n", + "Since $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 if the condition $-(n-1)+N \\ge m$ is satisfied. The condition can be written as $N \\ge n + m - 1$, which is satisfied!" + ] + }, + { + "cell_type": "markdown", + "id": "181cc549-9b2a-48ac-9077-71174cb49f1c", + "metadata": {}, + "source": [ + "# Option I: Bonus\n", + "\n", + "Let's revisit the linear convolution, shown in eq (2). Recall that the length of output is $n + m - 1$. The sliding dot product between $Q$ and $T$ is in $range(m-1, n)$ though. Therefore, if the goal is to calculate the sliding dot product, the computation can stop at index $n$. In other words, $idx$ can change from 0 to $n$ (exclusive). How about the N-Circular convolution in this case? If we set $N$ to $n$ instead of $n+m-1$, does the outcome have the correct answer for the slice in $range(m-1,n)$? We need to show that eq(2) and eq(8) give the same result for $m-1 \\le idx \\le n-1$ when $N=n$.\n", + "\n", + "**Case A: $j \\ge 0$**
\n", + "\n", + "Same as before!\n", + "\n", + "**Case B: $j < 0$**
\n", + "In this case, $Q'[j]=0$.\n", + "\n", + "In this case, $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$. Note that: \n", + "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge ((m-1) - (n-1)) + N$$ \n", + "\n", + "Since $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 if the condition $m-n + N \\ge m$ is satisfied. The condition can be written as $N \\ge n$, which is satisfied!" + ] + }, + { + "cell_type": "markdown", + "id": "03cbd90d-a7e4-466d-af15-81d9a1c4ce6c", + "metadata": {}, + "source": [ + "This shows that the sliding dot product of $Q$ and $T$ can be computed via $n$-Circular convolution rather than $(n+m-1)$-circular convolution" ] }, { From 520990368fc273da1ca82b8f4afdc1a50abd6e1b Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Thu, 3 Sep 2026 22:32:28 -0400 Subject: [PATCH 3/5] minor changes --- docs/overlap_add_math.ipynb | 62 ++++++++++++++++++++----------------- 1 file changed, 33 insertions(+), 29 deletions(-) diff --git a/docs/overlap_add_math.ipynb b/docs/overlap_add_math.ipynb index 8d0b9e7..fa43df6 100644 --- a/docs/overlap_add_math.ipynb +++ b/docs/overlap_add_math.ipynb @@ -5,7 +5,7 @@ "id": "e984e902-bf76-43ae-a5d5-55e59357b6bf", "metadata": {}, "source": [ - "# From Linear Convolution to Overlap-Add Method!" + "# From Linear Convolution to Overlap-Add Method" ] }, { @@ -13,7 +13,7 @@ "id": "e9725572-a39b-4fb5-992c-5f3d22a80797", "metadata": {}, "source": [ - "This notebook explains the math behind the overlap-add method in convolution. The content is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/)" + "This notebook explains the math behind linear convolution, circular convolution, and the overlap-add method in convolution. The content is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/)" ] }, { @@ -21,7 +21,7 @@ "id": "03f9d2f8-1560-4a5d-93a2-ec82fb422f6f", "metadata": {}, "source": [ - "Most textbooks use $x$ and $h$ notation in the convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. We would also like to introduce $Q'$, which is simply the reverse of $Q$." + "Most textbooks use $x$ and $h$ notation as the two signals in convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. Also, we would like to introduce $Q'$, which is simply the reverse of $Q$." ] }, { @@ -37,7 +37,7 @@ "id": "9ee474bd-1229-4749-80fe-f90c4a135931", "metadata": {}, "source": [ - "Consider two finte-length sequences $Q'$ (of length $m$) and $T$ (of length $n$). Their linear convolution, $C$, can be computed as follows:" + "Consider two finite-length sequences, $Q'$ and $T$, of lengths $m$ and $n$ respectively, where $m < n$. Their linear convolution, $C$, can be computed as follows:" ] }, { @@ -45,7 +45,7 @@ "id": "ad4a63b5-2d91-4875-b4c2-bc7ce11b6af9", "metadata": {}, "source": [ - "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$" + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]}, \\quad 0 \\le idx < n + m - 1$$" ] }, { @@ -55,7 +55,7 @@ "source": [ "where $T[.]$ and $Q'[.]$ are zeros for any index that is outside of their range. So, $T[i]$ is zero if $i$ is outside of the range $0 \\le i \\le n-1$. Hence, the equation above is equivalent to:\n", "\n", - "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]} $$" + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}, \\quad 0 \\le idx < n + m - 1$$" ] }, { @@ -63,7 +63,7 @@ "id": "593a8276-49cb-43d2-9300-8d29e67d3104", "metadata": {}, "source": [ - "Note that $Q'[idx-i]$ is not zero when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. \n", + "Note that $Q'[idx-i]$ can be a non-zero value when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. \n", "\n", "Therefore:\n", "\n", @@ -76,9 +76,7 @@ "id": "e638d3d7-cf2a-4355-aff3-a9192400defb", "metadata": {}, "source": [ - "This shows that the linear convolution $C$ has values at indices $\\set{0, 1, ..., n + m - 2}$, and it is 0 otherwise. Therefore, the length of output in linear convolution is $n+m-1$. \n", - "\n", - "Furthermore, let's check out the values for different indices:\n", + "This shows that the linear convolution $C$ has values at indices $\\set{0, 1, ..., n + m - 2}$, and it is 0 otherwise. Furthermore, let's check out the values at different indices:\n", "* $idx=0 \\implies C[0]=T[0]Q'[0]$\n", "* $idx=1 \\implies C[0]=T[0]Q'[1] + T[1]Q'[0]$\n", "* ...\n", @@ -89,15 +87,9 @@ "* ...\n", "* $idx=n+m-2 \\implies C[n+m-2]=T[n-1]Q'[(n+m-2)-(n-1)] = T[n-1]Q'[m-1]$\n", "\n", - "As observed, when $m-1 \\le idx \\le n-1$, $C[idx]$ becomes the dot product between a subsequnce of $T$ and $Q$, which is the reverse of $Q'$. Therefore, the `range(m-1,n)` gives the sliding dot product between $Q$ and $T$." - ] - }, - { - "cell_type": "markdown", - "id": "48c7f9d8-1b7c-4abc-886b-27a015b78ed6", - "metadata": {}, - "source": [ - "**How can we leverage this to speed up the computation of sliding dot product?**" + "As observed, when $m-1 \\le idx \\le n-1$, $C[idx]$ becomes the dot product between $Q$ (the reverse of $Q'$) and a subsequnce of $T$. Therefore, $range(m-1,n)$ gives the sliding dot product between $Q$ and $T$. \n", + "\n", + "The time complexity of calcularing the linear convolution is $O(nm)$. Can we do better? Yes!" ] }, { @@ -159,11 +151,11 @@ "id": "5ce64782-66d0-413b-a116-9ac3555966b6", "metadata": {}, "source": [ - "Let's break down the circular convolution equation:\n", + "Let's prove this! We start with breaking down the circular convolution equation:\n", "\n", - "$$C_{N}[idx] = \\sum_{i=0}^{m-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=m}^{n-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=n}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", + "$$C_{N}[idx] = \\sum_{i=0}^{n-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=n}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", "\n", - "$$C_{N}[idx] = \\sum_{i=0}^{m-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=m}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", + "In the second term, $i \\ge n$. Therefore, $\\tilde{T}[i]$ is zero. Hence, the equation above becomes:\n", "\n", "$$C_{N}[idx] = \\sum_{i=0}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$" ] @@ -173,9 +165,9 @@ "id": "a7e66371-4cd4-4673-957c-b7c548dcbfde", "metadata": {}, "source": [ - "Let's compare the linear convolution in eq(2) and the N-circular convolution in eq(9). To show that they give the same value, we need to compare $Q'[idx-i]$ with $\\tilde{Q'}[(idx-i)_{N}]$.\n", + "Let's compare the linear convolution in eq(2) and the N-circular convolution in eq(7). To show that they give the same value, we need to compare $Q'[idx-i]$ with $\\tilde{Q'}[(idx-i)_{N}]$.\n", "\n", - "Let $j$ denote $idx-i$. So, all we need to do is to compare $Q'[j]$ with $\\tilde{Q'}[(j)_{N}]$.\n", + "Let's define: $j = idx-i$. So, all we need to do is to compare $Q'[j]$ with $\\tilde{Q'}[(j)_{N}]$.\n", "\n", "**Case A: $j \\ge 0$**
\n", "\n", @@ -197,9 +189,9 @@ "id": "181cc549-9b2a-48ac-9077-71174cb49f1c", "metadata": {}, "source": [ - "# Option I: Bonus\n", + "# Option I: Enhanced\n", "\n", - "Let's revisit the linear convolution, shown in eq (2). Recall that the length of output is $n + m - 1$. The sliding dot product between $Q$ and $T$ is in $range(m-1, n)$ though. Therefore, if the goal is to calculate the sliding dot product, the computation can stop at index $n$. In other words, $idx$ can change from 0 to $n$ (exclusive). How about the N-Circular convolution in this case? If we set $N$ to $n$ instead of $n+m-1$, does the outcome have the correct answer for the slice in $range(m-1,n)$? We need to show that eq(2) and eq(8) give the same result for $m-1 \\le idx \\le n-1$ when $N=n$.\n", + "Let's revisit the linear convolution, shown in eq (2). Recall that the length of output is $n + m - 1$. The sliding dot product between $Q$ and $T$ is in $range(m-1, n)$ though. Therefore, if the goal is to calculate the sliding dot product, the computation can stop at index $n$. In other words, $idx$ can change from 0 to $n$ (exclusive). How about the N-Circular convolution in this case? If we set $N$ to $n$ instead of $n+m-1$, does the outcome have the correct answer for the slice in $range(m-1,n)$? We need to show that eq(2) and eq(7) give the same result for $m-1 \\le idx \\le n-1$ when $N=n$.\n", "\n", "**Case A: $j \\ge 0$**
\n", "\n", @@ -209,9 +201,13 @@ "In this case, $Q'[j]=0$.\n", "\n", "In this case, $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$. Note that: \n", - "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge ((m-1) - (n-1)) + N$$ \n", + "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge ((m-1) - (n-1)) + N \\ge m - n + N$$ \n", "\n", - "Since $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 if the condition $m-n + N \\ge m$ is satisfied. The condition can be written as $N \\ge n$, which is satisfied!" + "Recall that $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 as long as $j+N \\ge m$. the condition $m-n + N \\ge m$ is satisfied. Earlier, we showed that $j + N \\ge m - n + N$. Therefore, we just need to make sure $m - n + N$ is at least $m$. \n", + "\n", + "$$m - n + N \\ge m \\implies N \\ge n$$ \n", + "\n", + "And $N \\ge n$ is satisfied!" ] }, { @@ -219,7 +215,15 @@ "id": "03cbd90d-a7e4-466d-af15-81d9a1c4ce6c", "metadata": {}, "source": [ - "This shows that the sliding dot product of $Q$ and $T$ can be computed via $n$-Circular convolution rather than $(n+m-1)$-circular convolution" + "We just showed that the sliding dot product of $Q$ and $T$ can be computed via $n$-Circular convolution. Computing this in the frequency domain will have the time complexity of $O(nlogn)$." + ] + }, + { + "cell_type": "markdown", + "id": "d21fceee-2bab-4d01-9bee-83f73fe53e08", + "metadata": {}, + "source": [ + "Can we do better? See the next approach!" ] }, { From 24ab0d0a76a0e8e483d06a9188942a0de2744dfa Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Mon, 7 Sep 2026 23:16:19 -0400 Subject: [PATCH 4/5] revise objective section and structure --- docs/overlap_add_math.ipynb | 307 ++++++++++++++++++++++++++---------- 1 file changed, 220 insertions(+), 87 deletions(-) diff --git a/docs/overlap_add_math.ipynb b/docs/overlap_add_math.ipynb index fa43df6..ecc5cf0 100644 --- a/docs/overlap_add_math.ipynb +++ b/docs/overlap_add_math.ipynb @@ -2,18 +2,23 @@ "cells": [ { "cell_type": "markdown", - "id": "e984e902-bf76-43ae-a5d5-55e59357b6bf", + "id": "5c9355b9-0bd4-4740-9721-7394a19f4588", "metadata": {}, "source": [ - "# From Linear Convolution to Overlap-Add Method" + "# Objective" ] }, { "cell_type": "markdown", - "id": "e9725572-a39b-4fb5-992c-5f3d22a80797", + "id": "199eb655-fde9-4ee1-82e4-55bad5892f22", "metadata": {}, "source": [ - "This notebook explains the math behind linear convolution, circular convolution, and the overlap-add method in convolution. The content is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/)" + "The paper [\"Matrix Profile I\"](https://www.cs.ucr.edu/~eamonn/PID4481997_extend_Matrix%20Profile_I.pdf) uses the [MASS algorithm](https://www.cs.unm.edu/~mueen/FastestSimilaritySearch.html) to compute the distance between a query $Q$ and every subsequence of length $len(Q)$ in $T$. As part of this algorithm, the sliding dot product (sdp) is calculated. The [MASS paper](https://link.springer.com/epdf/10.1007/s10618-024-01005-2?sharing_token=067pAxnaLDvz89q1n_GJt_e4RwlQNchNByi7wbcMAY6HNwOWuMxQSNE3HcKcuL8siHB8L8krJpchQVaGvicUgoegxxV7BWgaU4Y9evg1FU1LGtxlvM9A5UrxrtkqDHyvoOPk_ttVNps-6-LSRf1vuRgFWJI7qrkqWGitjsXfRsw%3D) provided different versions of MASS for computing the sliding dot product. At their core, they all come down to computing one/more convolution(s) via the frequency domain. To understand how convolution comes into the picture and how it is used for computing the sliding dot product, this notebook is created to answer the following questions:\n", + "\n", + "1. **How is the sliding dot product (sdp) related to (linear) convolution?** This can help us discover the relationship between sdp and linear convolution. Convolution is a well-studied area and understanding this relationship allows us to leverage convolution methods to compute the sliding dot product.\n", + "2. **How can we compute the convolution faster?** Once we understand the relationship between sdp and convolution, we can focus on learning/using method for faster convolution.\n", + "3. **Can we reduce the output size, and hence the computational load without losing the sdp?** Once we know how to get faster convolution, we can try to enhance its performance further by avoiding unnecessary calculation.\n", + "4. **What to do if $T$ is very long? Overlap-add method!** This is where we learn a divide-and-conquer algorithm on convolution, a useful method when input $T$ is large" ] }, { @@ -21,7 +26,9 @@ "id": "03f9d2f8-1560-4a5d-93a2-ec82fb422f6f", "metadata": {}, "source": [ - "Most textbooks use $x$ and $h$ notation as the two signals in convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. Also, we would like to introduce $Q'$, which is simply the reverse of $Q$." + "The content of this notebook is based on the book \"Discrete-Time Signal Processing\" by Alan V. Oppenheim and Ronald W. Schafer. The book is [publicly available in MIT OpenCourseWare](https://ocw.mit.edu/courses/res-6-dtsp-discrete-time-signal-processing/resources/mitres_6-dtsp_s26_thirdedition_pdf/). Most textbooks use $x$ and $h$ notation as the two signals in convolution. However, in this notebook, we would like to use $Q$ and $T$ as those are the name of arrays in the realm of STUMPY. Also, we would like to introduce $Q'$, which is simply the reverse of $Q$. So, if $Q$ has length $m$ (indices $0,\\ldots,m-1$), then its reverse $Q'$ is defined by\n", + "\n", + "$$Q'[j] = Q[m-1-j], \\quad 0 \\le j \\le m-1,$$" ] }, { @@ -29,7 +36,15 @@ "id": "3dda604a-3a3e-42f4-90c7-7b01ada16337", "metadata": {}, "source": [ - "# Linear Convolution of Two Finite-Length Sequences" + "# 1. How is the sliding dot product(sdp) related to (linear) convolution?" + ] + }, + { + "cell_type": "markdown", + "id": "877af78c-0026-44c0-8295-0f75db075d97", + "metadata": {}, + "source": [ + "### 1.1 Linear Convolution of Two Finite-Length Sequences $Q'$ and $T$" ] }, { @@ -37,7 +52,7 @@ "id": "9ee474bd-1229-4749-80fe-f90c4a135931", "metadata": {}, "source": [ - "Consider two finite-length sequences, $Q'$ and $T$, of lengths $m$ and $n$ respectively, where $m < n$. Their linear convolution, $C$, can be computed as follows:" + "Consider two finite-length sequences, $Q'$ and $T$, of lengths $m$ and $n$, respectively, where $m < n$. Their linear convolution, $C$, can be computed as follows:" ] }, { @@ -45,7 +60,7 @@ "id": "ad4a63b5-2d91-4875-b4c2-bc7ce11b6af9", "metadata": {}, "source": [ - "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]}, \\quad 0 \\le idx < n + m - 1$$" + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]}$$" ] }, { @@ -53,9 +68,11 @@ "id": "8028ee30-3244-4c79-86ca-22e2b8876694", "metadata": {}, "source": [ - "where $T[.]$ and $Q'[.]$ are zeros for any index that is outside of their range. So, $T[i]$ is zero if $i$ is outside of the range $0 \\le i \\le n-1$. Hence, the equation above is equivalent to:\n", + "where $T[.]$ and $Q'[.]$ are zeros for any index that is outside of their range. The linear convolution can be seen as a flip and slide operation, and at each index $idx$, the output is the sum of element-wise product of overlapping elements. So, as long as there is at least one overlapping element, there can be an output. For instance, at index $idx=0$, there is only one overlapping element between $T$ and $Q'$, and that $T[0]$ and $Q'[0]$. And, $C[0] = T[0]Q'[0]$. Therefore, **the linear convolution is NOT the same as the sliding dot product** since sliding dot product is for cases when the overlap covers the full query and not just one element. \n", + "\n", + "As stated earlier, $T[i]$ is zero if $i$ is outside of the range $0 \\le i \\le n-1$. Hence, the boundaries of the sum operation in the previous equation can be revised as follows:\n", "\n", - "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}, \\quad 0 \\le idx < n + m - 1$$" + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}$$" ] }, { @@ -63,12 +80,23 @@ "id": "593a8276-49cb-43d2-9300-8d29e67d3104", "metadata": {}, "source": [ - "Note that $Q'[idx-i]$ can be a non-zero value when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. \n", + "Note that $Q'[idx-i]$ can be a non-zero value when $0 \\le idx-i \\le m-1$, or equivalently $i \\le idx \\le i+m-1$. So, as the index $i$ changes from $0$ to $n-1$, the interval $i \\le idx \\le i+m-1$ changes accordingly, and the range of $idx$ becomes the union of these intervals, i.e.\n", "\n", - "Therefore:\n", + "$$\\bigcup_{i=0}^{n-1} [\\,i,\\ i+m-1\\,]$$\n", + "\n", + "So:\n", + "* min(idx) is coming from the lower bound of interval $[\\,i,\\ i+m-1\\,]$ when $i=0$. This gives: $min(idx)=0$\n", + "* max (idx) is coming from the upper bound of interval $[\\,i,\\ i+m-1\\,]$ when $i=n-1$. This gives: $max(idx) = (n-1)+m-1 = n+m-2$\n", "\n", - "* $idx \\ge i, \\quad i \\ge 0 \\implies idx \\ge 0$\n", - "* $idx \\le i+m-1, \\quad i \\le n-1 \\implies idx \\le n+m-2$" + "So, when $idx$ changes from $0$ to $n+m-2$, there exist an $i$ such that both $T[i]$ and $Q'[idx-i]$ can be non-zero values. This shows that the linear convolution $C$ can have values at indices $\\set{0, 1, ..., n + m - 2}$. " + ] + }, + { + "cell_type": "markdown", + "id": "ecb14f7d-cde6-42db-aece-2fb40536d68c", + "metadata": {}, + "source": [ + "### 1.2 Sliding dot product (of $Q$ and $T$) as a slice of linear convolution (of $Q'$ and $T$)" ] }, { @@ -76,20 +104,33 @@ "id": "e638d3d7-cf2a-4355-aff3-a9192400defb", "metadata": {}, "source": [ - "This shows that the linear convolution $C$ has values at indices $\\set{0, 1, ..., n + m - 2}$, and it is 0 otherwise. Furthermore, let's check out the values at different indices:\n", + "Let's check out the values of linear convolution at different indices. For convenience, the linear convolution equation is repeated below:\n", + "\n", + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}$$\n", + "\n", + "> **NOTE** Let $T_{i}$ denote the sequence $\\{T[i], T[i+1], ..., T[i+m-1]\\}$\n", + "\n", "* $idx=0 \\implies C[0]=T[0]Q'[0]$\n", "* $idx=1 \\implies C[0]=T[0]Q'[1] + T[1]Q'[0]$\n", "* ...\n", - "* $idx=m-1 \\implies C[m-1]=T[0]Q'[m-1] + T[1]Q'[m-2] + ... + T[m-1]Q'[0] = T_{0}.Q$\n", - "* $idx=m \\implies C[m]=T[1]Q'[m-1] + T[2]Q'[m-2] + ... + T[m]Q'[0] = T_{1}.Q$\n", + "* $\\textcolor{red}{idx=m-1 \\implies C[m-1]=T[0]Q'[m-1] + T[1]Q'[m-2] + ... + T[m-1]Q'[0] = T_{0}.Q}$\n", + "* $\\textcolor{red}{idx=m \\implies C[m]=T[1]Q'[m-1] + T[2]Q'[m-2] + ... + T[m]Q'[0] = T_{1}.Q}$\n", "* ...\n", - "* $idx=n-1 \\implies C[n-1]=T[n-m]Q'[(n-1)-(n-m)] + T[n-m+1]Q'[(n-1)-(n-m+1)] + ... + T[n-1]Q'[0] = T_{n-m+1}.Q$\n", + "* $\\textcolor{red}{idx=n-1 \\implies C[n-1]=T[n-m]Q'[(n-1)-(n-m)] + T[n-m+1]Q'[(n-1)-(n-m+1)] + ... + T[n-1]Q'[0] = T_{n-m}.Q}$\n", "* ...\n", "* $idx=n+m-2 \\implies C[n+m-2]=T[n-1]Q'[(n+m-2)-(n-1)] = T[n-1]Q'[m-1]$\n", "\n", - "As observed, when $m-1 \\le idx \\le n-1$, $C[idx]$ becomes the dot product between $Q$ (the reverse of $Q'$) and a subsequnce of $T$. Therefore, $range(m-1,n)$ gives the sliding dot product between $Q$ and $T$. \n", + "When $m-1 \\le idx \\le n-1$, $C[idx]$ actually becomes the dot product between $Q$ (the reverse of $Q'$) and $T_{idx-(m-1)}$. Therefore, the linear convolution of $Q'$ and $T$ contains the sliding dot product between $Q$ and $T$ in the slice in $range(m-1,n)$. \n", "\n", - "The time complexity of calcularing the linear convolution is $O(nm)$. Can we do better? Yes!" + "So, once we calculate the linear convolution between $Q'$ and $T$, we have the sliding dot product between $Q$ and $T$ for free! The time complexity of calculating the linear convolution is $O(nm)$. Can we do better? " + ] + }, + { + "cell_type": "markdown", + "id": "8475be6a-9515-43ae-a8d7-083b8db3af51", + "metadata": {}, + "source": [ + "# 2. How can we compute the convolution faster than O(nm)?" ] }, { @@ -97,7 +138,7 @@ "id": "52bcf566-3ac0-483d-81a3-57fa3b5766c3", "metadata": {}, "source": [ - "# Option I: Convert to Circular Convolution and use FFT-IFFT" + "### 2.1 Convert linear convolution to circular convolution (Step I)" ] }, { @@ -105,11 +146,11 @@ "id": "2c4176e2-60f3-410f-ad2d-6c1d29cca6cf", "metadata": {}, "source": [ - "Circular convolution is defined between two sequences that are both periodic and their period are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n", + "Circular convolution is defined between two sequences that are both periodic and their periods are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n", "\n", "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", "\n", - "where, $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$. $C_{N}$ represents one period of N-Circular convolution. This can also be computed via FFT, i.e. $C_{N} = IFFT( FFT(\\tilde{T}_{N}) \\times FFT(\\tilde{Q'}_{N}) )$." + "where, $C_{N}$ represents one period of convolution output. $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$, and the index $(idx-i)_{N}$ simply means \"$idx-i$ modulo $N$\", i.e. $(idx-i) \\% N$, and it is always between $0$ and $N-1$." ] }, { @@ -117,45 +158,50 @@ "id": "e10c5304-a1a1-4cbd-9b26-266e229bdc57", "metadata": {}, "source": [ - "If there is a way to compute the linear convolution via circular convolution, then we can use FFT to compute the linear convolution. Actually, there is a way! The linear convolution between $Q'$ and $T$ is equivalent to the N-circular convolution between $\\tilde{Q'}$ and $\\tilde{T}$, where:\n", + "The time complexity of calculating $C_{N}$ using the equation above is $O(N^{2})$. However, as will be shown later in the next section, the time complexity can be reduced to $O(NlogN)$. However, this means nothing for sliding dot product unless we understand how the linear convolution can be computed via circular convolution and how $N$ is related to $n$ and/or $m$. \n", + "\n", + "We first prove that the linear convolution between $Q'$ and $T$ is equivalent to $C_{N}$ between $\\tilde{Q'}$ and $\\tilde{T}$ if:\n", "\n", - "* $N \\ge n + m - 1$, to cover the $n+m-1$ data points of linear convolution\n", - "* $\\tilde{Q'}_{N}$ is $Q'$ but zero-padded with $N-m$ zeros\n", - "* $\\tilde{T}_{N}$ is $T$ but zero-padded with $N-n$ zeros" + "* $N = n + m - 1$\n", + "* $\\tilde{Q'}_{N}$, i.e. one period of $\\tilde{Q'}$, is $Q'$ but with $N-m$ zero padding\n", + "* $\\tilde{T}_{N}$, i.e. one period of $\\tilde{T}$, is $T$ but with $N-n$ zero padding" ] }, { "cell_type": "markdown", - "id": "b6529137-fcdc-4d65-aa93-da3746739f94", + "id": "5ce64782-66d0-413b-a116-9ac3555966b6", "metadata": {}, "source": [ - "So:\n", + "**Proof**\n", "\n", - "$$\n", + "> **Note:**
\n", + "> Keep in mind that $N = n + m - 1$, and:\n", + "> \n", + "> $$\n", "\\tilde{Q'}_{N}[i] = \\begin{cases} \n", " Q'[i] & \\text{if } 0 \\le i \\le m-1 \\\\\n", " 0 & \\text{if } m \\le i \\le N-1\n", " \\end{cases}\n", "$$\n", - "\n", - "$$\n", + "> \n", + "> $$\n", "\\tilde{T}_{N}[i] = \\begin{cases} \n", " T[i] & \\text{if } 0 \\le i \\le n-1 \\\\\n", " 0 & \\text{if } n \\le i \\le N-1\n", " \\end{cases}\n", - "$$" - ] - }, - { - "cell_type": "markdown", - "id": "5ce64782-66d0-413b-a116-9ac3555966b6", - "metadata": {}, - "source": [ - "Let's prove this! We start with breaking down the circular convolution equation:\n", + "$$\n", + "\n", + "\n", + "\n", + "We start with rewriting the circular convolution equation by breaking the summation into two parts: \n", "\n", "$$C_{N}[idx] = \\sum_{i=0}^{n-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=n}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", "\n", - "In the second term, $i \\ge n$. Therefore, $\\tilde{T}[i]$ is zero. Hence, the equation above becomes:\n", + "Note that:\n", + "* In the first term, $0 \\le i \\le n-1$, and hence $\\tilde{T}[i]$ is the same as $T[i]$.\n", + "* In the second term, $n \\le i \\le N-1$, and hence $\\tilde{T}[i]$ is zero.\n", + "\n", + "Hence, the equation becomes:\n", "\n", "$$C_{N}[idx] = \\sum_{i=0}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$" ] @@ -165,65 +211,155 @@ "id": "a7e66371-4cd4-4673-957c-b7c548dcbfde", "metadata": {}, "source": [ - "Let's compare the linear convolution in eq(2) and the N-circular convolution in eq(7). To show that they give the same value, we need to compare $Q'[idx-i]$ with $\\tilde{Q'}[(idx-i)_{N}]$.\n", + "Recall that the linear convolution equation is:\n", + "\n", + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}$$\n", + "\n", + "So, to show that $C_{N}$, from circular convolution, gives the same value as the linear convolution, we only need to show that $\\tilde{Q'}[(idx-i)_{N}]$ and $Q'[idx-i]$ give the same value when $N=n+m-1$. Let $j$ denote $idx-i$. So, all we need to do is to prove that $\\tilde{Q'}[(j)_{N}]$ and $Q'[j]$ give the same value for different values of $j$.\n", + "\n", + "**Case I: $j \\ge 0$**
\n", + "\n", + "Let's compute the upper bound for $j$. Recall that: \n", + "* $j=idx-i$\n", + "* $0 \\le idx \\le N-1$ (according to range of indices of one period in circular convolution)\n", + "* $0 \\le i \\le n-1$ (according to range of indices of elements in $T$)\n", + "\n", + "Therefore,\n", "\n", - "Let's define: $j = idx-i$. So, all we need to do is to compare $Q'[j]$ with $\\tilde{Q'}[(j)_{N}]$.\n", + "$$max(j) = max(idx-i) = max(idx) + max(-i) = max(idx) - min(i) = (N-1) - 0 = N-1$$ \n", "\n", - "**Case A: $j \\ge 0$**
\n", + "So, $j$ changes from 0 to $N-1$. hence, $j$ modulo $N$ is still $j$. Therefore: $\\tilde{Q'}[(j)_{N}]=Q'[j]$. Proof is now complete for this case.\n", "\n", - "In this case, $max(j) = max(idx) - min(i) = (N-1) - 0 = N-1 \\implies \\tilde{Q'}[(j)_{N}]=Q'_{N}[j]$\n", + "**Case II: $j < 0$**
\n", "\n", + "In this case, $Q'[j]=0$. So, all we need to do is to prove that $\\tilde{Q'}[(j)_{N}]$ is zero as well. \n", "\n", - "**Case B: $j < 0$**
\n", + "Let's compute the lower bound for $j$. Recall that: \n", + "* $j=idx-i$\n", + "* $0 \\le idx \\le N-1$\n", + "* $0 \\le i \\le n-1$\n", "\n", - "In this case, $Q'[j]=0$.\n", + "Therefore,\n", "\n", - "In this case, $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$. Note that: \n", - "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge (0 - (n-1)) + N$$ \n", + "$$min(j) = min(idx-i) = min(idx) + min(-i) = min(idx) - max(i) = 0 - (n-1)$$\n", "\n", - "Since $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 if the condition $-(n-1)+N \\ge m$ is satisfied. The condition can be written as $N \\ge n + m - 1$, which is satisfied!" + "So, in this case, $-(n-1) \\le j \\le -1$. \n", + "\n", + "Since $\\tilde{Q'}[(j)_{N}] == \\tilde{Q'}[(j+N)_{N}]$, we can compute the latter instead of the former. To do that, we first need to figure out the range for $j+N$.\n", + "\n", + "$$ -(n-1) \\le j \\le -1 \\implies -(n-1)+N \\le j+N \\le -1+N $$\n", + "\n", + "Let's plug the value of $N=n+m-1$ into $-(n-1)+N$:\n", + "\n", + "$$(n+m-1)-(n-1) \\le j+N \\le N-1 \\implies m \\le j+N \\le N-1$$\n", + "\n", + "Since $j+N$ is between $m$ and $N-1$, its value modulo $N$ becomes $j+N$ again, meaning $\\tilde{Q'}[(j)_{N}] == \\tilde{Q'}[(j+N)_{N}] = \\tilde{Q'}[j+N]$. Since the index $j+N$ is $\\ge m$, the value $\\tilde{Q'}[j+N]$ becomes 0 as the element resides in the zero-padding part. Proof is now complete for this case." ] }, { "cell_type": "markdown", - "id": "181cc549-9b2a-48ac-9077-71174cb49f1c", + "id": "66a707c1-8004-4fe6-b954-ffe039eda441", "metadata": {}, "source": [ - "# Option I: Enhanced\n", - "\n", - "Let's revisit the linear convolution, shown in eq (2). Recall that the length of output is $n + m - 1$. The sliding dot product between $Q$ and $T$ is in $range(m-1, n)$ though. Therefore, if the goal is to calculate the sliding dot product, the computation can stop at index $n$. In other words, $idx$ can change from 0 to $n$ (exclusive). How about the N-Circular convolution in this case? If we set $N$ to $n$ instead of $n+m-1$, does the outcome have the correct answer for the slice in $range(m-1,n)$? We need to show that eq(2) and eq(7) give the same result for $m-1 \\le idx \\le n-1$ when $N=n$.\n", - "\n", - "**Case A: $j \\ge 0$**
\n", - "\n", - "Same as before!\n", - "\n", - "**Case B: $j < 0$**
\n", - "In this case, $Q'[j]=0$.\n", + "We just proved that the linear convolution can be computed in the form of circular convolution when the arrays have certain lengths and are zero-padded properly." + ] + }, + { + "cell_type": "markdown", + "id": "8648ff47-4d67-4bec-a91e-4ec6d799cc0e", + "metadata": {}, + "source": [ + "### 2.2 Compute Circular Convolution in the Frequency Domain (Step II)" + ] + }, + { + "cell_type": "markdown", + "id": "a6b6ff9c-d1e2-4422-ab1b-e7f4915ea0f9", + "metadata": {}, + "source": [ + "The previous part showed that linear convolution between $Q'$ and $T$ can be computed via circular convolution when both arrays are zero-padded till they reach the length $N=n+m-1$, where $n$ and $m$ are the length of $T$ and $Q'$, respectively. However, the computation via the formula provided in the previous section is not faster. In fact, it has the time complexity of $O(N^2)$! So, why did we go through all that trouble to show that linear convolution can be computed via circular convolution? This is because the circular convolution can be computed more efficiently when it is computed with the help of Fast Fourier Transform. So, instead of computing the circular convolution in time domain, i.e.\n", "\n", - "In this case, $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$. Note that: \n", - "$$j+N \\ge min(j)+N \\ge (min(idx)-max(i)) + N \\ge ((m-1) - (n-1)) + N \\ge m - n + N$$ \n", + "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", "\n", - "Recall that $Q'_{N}$ is padded with zeros for indices that are at least $m$, therefore the value $\\tilde{Q'}[(j)_{N}] = Q'_{N}[j+N]$ can be 0 as long as $j+N \\ge m$. the condition $m-n + N \\ge m$ is satisfied. Earlier, we showed that $j + N \\ge m - n + N$. Therefore, we just need to make sure $m - n + N$ is at least $m$. \n", + "we can compute it by taking it to the frequency domain:\n", "\n", - "$$m - n + N \\ge m \\implies N \\ge n$$ \n", + "$$C_{N} = IFFT\\left(\n", + "FFT(\\tilde{T}_{N})\n", + "\\cdot\n", + "FFT(\\tilde{Q'}_{N})\n", + "\\right)$$\n", "\n", - "And $N \\ge n$ is satisfied!" + "and its time complexity becomes $O(NlogN)$." ] }, { "cell_type": "markdown", - "id": "03cbd90d-a7e4-466d-af15-81d9a1c4ce6c", + "id": "64b8a7ee-f1dd-4494-a9a4-4ba3293a3912", "metadata": {}, "source": [ - "We just showed that the sliding dot product of $Q$ and $T$ can be computed via $n$-Circular convolution. Computing this in the frequency domain will have the time complexity of $O(nlogn)$." + "So, if I can compute the circular convolution faster, it means that I can calculate the linear convolution faster, and therefore I can obtain the sliding dot product faster than before!" ] }, { "cell_type": "markdown", - "id": "d21fceee-2bab-4d01-9bee-83f73fe53e08", + "id": "181cc549-9b2a-48ac-9077-71174cb49f1c", "metadata": {}, "source": [ - "Can we do better? See the next approach!" + "# 3. Can we reduce the output size, and hence the computational load without losing the sdp?\n", + "\n", + "Let's revisit the linear convolution:\n", + "\n", + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}, \\quad 0 \\le idx \\le n+m-2$$\n", + "\n", + "\n", + "\n", + "Recall that the length of output is $n + m - 1$. The sliding dot product between $Q$ and $T$ is in $range(m-1, n)$ though. Therefore, if all we care about is the sliding dot product, the computation can stop at index $n$. In other words, $idx$ (of linear convolution) can be from $0$ to $n-1$. How about the Circular convolution? If we accordingly set $N$ to $n$ instead of $n+m-1$, does the slice in $range(m-1,n)$ still reflect the values of sdp? In other words, we need to show that\n", + "\n", + "$$ C[idx] = \\sum_{i=0}^{n-1}{T[i] \\times Q'[idx-i]}, \\quad m-1 \\le idx \\le n-1$$\n", + "\n", + "and\n", + "\n", + "$$C_{N=n}[idx] = \\sum_{i=0}^{n-1}T[i] \\times \\tilde{Q'}[(idx-i)_{n}], \\quad m-1 \\le idx \\le n-1$$\n", + "\n", + "give the same result. To prove this, we just need to show that $Q'[idx-i]$ and $\\tilde{Q'}[(idx-i)_{n}]$ have the same value when $m-1 \\le idx \\le n-1$.
\n", + "\n", + "\n", + "**Proof:** Let $j$ denote $idx-i$. We need to check two different cases:\n", + "\n", + "**Case I: $j \\ge 0$**
\n", + "\n", + "Let's compute the upper bound for $j$. Recall that: \n", + "\n", + "* $j=idx-i$\n", + "* $m-1 \\le idx \\le n-1$\n", + "* $0 \\le i \\le n-1$\n", + "\n", + "$$ j \\le max(j)=max(idx-i)=max(idx)+max(-i)=max(idx) - min(i)=(n-1)-0=n-1$$\n", + "\n", + "Since $0 \\le j \\le n-1$, $\\tilde{Q'}[(j)_{n}]$ is the same as $Q'[j]$. Proof is now complete for this case.\n", + "\n", + "**Case II: $j < 0$**
\n", + "Note that $Q'[j]=0$ for $j < 0$. So, we need to show $\\tilde{Q'}[(j)_{n}]$ becomes zero as well in this case.\n", + "\n", + "Let's compute the lower bound for $j$. Recall that:\n", + "\n", + "* $j=idx-i$\n", + "* $m-1 \\le idx \\le n-1$\n", + "* $0 \\le i \\le n-1$\n", + "\n", + " \n", + "* $\\tilde{Q'}[(j)_{n}] = \\tilde{Q'}[(j+n)_{n}]$. Note that: \n", + "$$j+n \\ge min(j)+n \\ge min(idx-i)+n \\ge min(idx)+min(-i) + n \\ge min(idx) - max(i) + n \\ge ((m-1) - (n-1)) + n \\ge m$$ \n", + "\n", + "So, when $j<0$, then: $m \\le j+n \\le n-1$. Therefore $\\tilde{Q'}[(j+n)_{n}]$ simply becomes $\\tilde{Q'}[j+n]$, and the value is zero as the index falls into the padded zeros. Therefore, $Q'[j]=\\tilde{Q'}[j+n]=0$. Proof is now complete for this case. " + ] + }, + { + "cell_type": "markdown", + "id": "03cbd90d-a7e4-466d-af15-81d9a1c4ce6c", + "metadata": {}, + "source": [ + "We just showed that the sliding dot product can be computed via circular convolution with period $N=n$. As noted in previous section, this can be computed in $O(nlogn)$ which is faster than our baseline $O(nm)$ unless $m$ is small." ] }, { @@ -231,8 +367,9 @@ "id": "97fcc328-35f9-4b66-9598-1b5d8f1f97ac", "metadata": {}, "source": [ - "# Option II: Overlap-add method\n", - "This is a divide and conquer algorithm. It applies divide on \"linear convolution\" and conquer each via circular convolution (See: Option I)" + "# 4. What to do if $T$ is a long sequence? Use Overlap-add method!\n", + "\n", + "It applies a divide-and-conquer algorithm on convolution. We first need to learn about the \"linearity\" property in linear convolution and how it allows us to divide the problem into similar problems but with smaller sizes. Then, we learn how we can use circular convolution on those smaller problems, and combine their outputs." ] }, { @@ -240,7 +377,7 @@ "id": "1e45ccf7-1564-4c41-b32e-4de9b1b201a3", "metadata": {}, "source": [ - "## Linearity in Linear Convolution" + "### 4.1 The \"Linearity\" property in Linear Convolution" ] }, { @@ -248,15 +385,15 @@ "id": "a084578e-cb20-437c-a70c-3e33f8da8ac3", "metadata": {}, "source": [ - "Suppose the array $T$ can be written as $T_{1} + T_{2}$. In other words: $T[i] = T_{1}[i] + T_{2}[i]$, then:\n", + "Suppose the array $T$ can be written as $T^{(1)} + T^{(2)}$. In other words: $T[i] = T^{(1)}[i] + T^{(2)}[i]$, then:\n", "\n", "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T[i] \\times Q'[idx-i]} $$\n", "\n", - "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{(T_{1}[i] \\times Q'[idx-i] + T_{2}[i] \\times Q'[idx-i])} $$\n", + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{(T^{(1)}[i] \\times Q'[idx-i] + T^{(2)}[i] \\times Q'[idx-i])} $$\n", "\n", - "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T_{1}[i] \\times Q'[idx-i]} + \\sum_{i=-\\infty}^{\\infty}{T_{2}[i] \\times Q'[idx-i]}$$\n", + "$$ C[idx] = \\sum_{i=-\\infty}^{\\infty}{T^{(1)}[i] \\times Q'[idx-i]} + \\sum_{i=-\\infty}^{\\infty}{T^{(2)}[i] \\times Q'[idx-i]}$$\n", "\n", - "$$ C[idx] = C_{1}[idx] + C_{2}[idx] $$" + "$$ C[idx] = C^{(1)}[idx] + C^{(2)}[idx] $$" ] }, { @@ -264,22 +401,18 @@ "id": "1f5be13b-82a7-46d9-9cb1-229c878e67f8", "metadata": {}, "source": [ - "This shows that a linear convolution has the property of \"linearity\". Now, we use this property to show that we can compute the linear convolution between $Q'$ and long $T$ by breaking $T$ into smaller parts." + "This shows that a linear convolution has the property of \"linearity\". Now, we use this property to show that we can compute the linear convolution between $Q'$ and $T$ by breaking $T$ into smaller parts." ] }, - { - "cell_type": "markdown", - "id": "fecee5fc-3c10-4c2e-8c66-bce17c3650bd", - "metadata": {}, - "source": [] - }, { "cell_type": "code", "execution_count": null, - "id": "acf5516d-e52d-457c-a8c5-f4d1b7ae8f1c", + "id": "fecee5fc-3c10-4c2e-8c66-bce17c3650bd", "metadata": {}, "outputs": [], - "source": [] + "source": [ + "#WIP" + ] } ], "metadata": { From 083c8d8f75e81f85904cac2f5c2cbd1c0a366ee8 Mon Sep 17 00:00:00 2001 From: NimaSarajpoor Date: Mon, 7 Sep 2026 23:23:04 -0400 Subject: [PATCH 5/5] minor changes --- docs/overlap_add_math.ipynb | 45 +++---------------------------------- 1 file changed, 3 insertions(+), 42 deletions(-) diff --git a/docs/overlap_add_math.ipynb b/docs/overlap_add_math.ipynb index ecc5cf0..731676e 100644 --- a/docs/overlap_add_math.ipynb +++ b/docs/overlap_add_math.ipynb @@ -146,11 +146,11 @@ "id": "2c4176e2-60f3-410f-ad2d-6c1d29cca6cf", "metadata": {}, "source": [ - "Circular convolution is defined between two sequences that are both periodic and their periods are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$, and it can be computed as follows:\n", + "Circular convolution is defined between two sequences that are both periodic and their periods are the same, say $N$. Their circular convolution is also a periodic sequence, with period $N$. The following equation shows how it can be computed for the period in $range(0, N)$.\n", "\n", "$$C_{N}[idx] = \\sum_{i=0}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}], \\quad 0 \\le idx \\le N-1$$\n", "\n", - "where, $C_{N}$ represents one period of convolution output. $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequence with period $N$, and the index $(idx-i)_{N}$ simply means \"$idx-i$ modulo $N$\", i.e. $(idx-i) \\% N$, and it is always between $0$ and $N-1$." + "where, $C_{N}$ represents one period of convolution output. $\\tilde{Q'}$ and $\\tilde{T}$ are both periodic sequences, each with period $N$. The index $(idx-i)_{N}$ simply means \"$idx-i$ modulo $N$\", i.e. $(idx-i) \\% N$, and it is always between $0$ and $N-1$." ] }, { @@ -158,7 +158,7 @@ "id": "e10c5304-a1a1-4cbd-9b26-266e229bdc57", "metadata": {}, "source": [ - "The time complexity of calculating $C_{N}$ using the equation above is $O(N^{2})$. However, as will be shown later in the next section, the time complexity can be reduced to $O(NlogN)$. However, this means nothing for sliding dot product unless we understand how the linear convolution can be computed via circular convolution and how $N$ is related to $n$ and/or $m$. \n", + "The time complexity of calculating $C_{N}$ using the equation above is $O(N^{2})$. However, as will be shown later in the next section, the time complexity can be reduced to $O(NlogN)$. However, this means nothing for sliding dot product unless we understand how circular convolution can be used to compute the linear convolution, and how $N$ is related to $n$ and/or $m$. \n", "\n", "We first prove that the linear convolution between $Q'$ and $T$ is equivalent to $C_{N}$ between $\\tilde{Q'}$ and $\\tilde{T}$ if:\n", "\n", @@ -167,45 +167,6 @@ "* $\\tilde{T}_{N}$, i.e. one period of $\\tilde{T}$, is $T$ but with $N-n$ zero padding" ] }, - { - "cell_type": "markdown", - "id": "5ce64782-66d0-413b-a116-9ac3555966b6", - "metadata": {}, - "source": [ - "**Proof**\n", - "\n", - "> **Note:**
\n", - "> Keep in mind that $N = n + m - 1$, and:\n", - "> \n", - "> $$\n", - "\\tilde{Q'}_{N}[i] = \\begin{cases} \n", - " Q'[i] & \\text{if } 0 \\le i \\le m-1 \\\\\n", - " 0 & \\text{if } m \\le i \\le N-1\n", - " \\end{cases}\n", - "$$\n", - "> \n", - "> $$\n", - "\\tilde{T}_{N}[i] = \\begin{cases} \n", - " T[i] & \\text{if } 0 \\le i \\le n-1 \\\\\n", - " 0 & \\text{if } n \\le i \\le N-1\n", - " \\end{cases}\n", - "$$\n", - "\n", - "\n", - "\n", - "We start with rewriting the circular convolution equation by breaking the summation into two parts: \n", - "\n", - "$$C_{N}[idx] = \\sum_{i=0}^{n-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}] + \\sum_{i=n}^{N-1} \\tilde{T}[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$\n", - "\n", - "Note that:\n", - "* In the first term, $0 \\le i \\le n-1$, and hence $\\tilde{T}[i]$ is the same as $T[i]$.\n", - "* In the second term, $n \\le i \\le N-1$, and hence $\\tilde{T}[i]$ is zero.\n", - "\n", - "Hence, the equation becomes:\n", - "\n", - "$$C_{N}[idx] = \\sum_{i=0}^{n-1} T[i] \\times \\tilde{Q'}[(idx-i)_{N}]$$" - ] - }, { "cell_type": "markdown", "id": "a7e66371-4cd4-4673-957c-b7c548dcbfde",