added dft play, renamed notebooks folder, since it now contains more …

…than just examples, and bumped version number because change of readme to link other places will merit a rerelease

README.md

@@ -28,7 +28,7 @@ You should now be able to
 >>> yth_n = np.sin(th_n) # must be periodic on the domain
 >>> dyth_n = fourier_deriv(yth_n, 1)
 ```
-For further usage examples, see the Jupyter notebooks in the examples folder: [Chebyshev](https://github.com/pavelkomarov/spectral-derivatives/blob/main/examples/chebyshev.ipynb) and [Fourier](https://github.com/pavelkomarov/spectral-derivatives/blob/main/examples/fourier.ipynb)
+For further usage examples, see the Jupyter notebooks in the notebooks folder: [Chebyshev](https://github.com/pavelkomarov/spectral-derivatives/blob/main/notebooks/chebyshev.ipynb) and [Fourier](https://github.com/pavelkomarov/spectral-derivatives/blob/main/notebooks/fourier.ipynb)
 
 ## References
 

notebooks/dft_play.ipynb

@@ -0,0 +1,208 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "738dfead-c1af-4a2e-ba1a-e523cb5190ce",
+   "metadata": {},
+   "source": [
+    "# Playing with the DFT\n",
+    "\n",
+    "The Discrete Fourier Transform definition can be twisted into several equivalent forms, which may be worth understanding. I use [very similar manipulations to prove $Y_{M-k} = Y_k$ in the math](https://pavelkomarov.com/spectral-derivatives/math.pdf)."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "18531c04-8603-4dc2-b026-63a2b1aef5ea",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "\n",
+    "def FFT(k):\n",
+    "    \"\"\"just to show DFT(k) definitely == FFT(k) == a_k\"\"\"\n",
+    "    return np.fft.fft(x)[k]\n",
+    "\n",
+    "def DFT(k):\n",
+    "\t\"\"\"standard way to find DFT coefficients\"\"\"\n",
+    "\ta_k = 0\n",
+    "\tfor n in range(len(x)):\n",
+    "\t\ta_k += x[n] * np.exp(-2*np.pi/len(x) * 1j * n * k)\n",
+    "\treturn a_k\n",
+    "\n",
+    "def DFT_(k):\n",
+    "\t\"\"\"finding DFT coefficients by summing in the opposite direction, using -n\"\"\"\n",
+    "\ta_k = 0\n",
+    "\tfor n in range(1, len(x)+1):\n",
+    "\t\ta_k += x[-n] * np.exp(2*np.pi/len(x) * 1j * n * k)\n",
+    "\treturn a_k\n",
+    "\n",
+    "def DFT__(k, L=3):\n",
+    "\t\"\"\"finding DFT coefficients with arbitrary roll-around, a mix of the above two formulae\"\"\"\n",
+    "\ta_k = 0\n",
+    "\tfor n in range(len(x)-L+1, len(x)+1): # from M-L to M, covers first half\n",
+    "\t\ta_k += x[-n] * np.exp(2*np.pi/len(x) * 1j * n * k)\n",
+    "\tfor n in range(L, len(x)):\n",
+    "\t\ta_k += x[n] * np.exp(-2*np.pi/len(x) * 1j * n * k)\n",
+    "\treturn a_k\n",
+    "\n",
+    "def DFT_even(k):\n",
+    "\t\"\"\"a formula which equals the other three for *even* signals\"\"\"\n",
+    "\ta_k = 0\n",
+    "\tfor n in range(len(x)):\n",
+    "\t\ta_k += x[n] * np.exp(2*np.pi/len(x) * 1j * n * k)\n",
+    "\treturn a_k"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "d8478cfa-acbc-4573-8e40-4f09e25d8861",
+   "metadata": {},
+   "source": [
+    "Let's try it on an even signal. This makes the input vector palindromic: $[x_0, x_1, ... x_N, x_{-(N-1)}, ...x_{-1}]$"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "0aaab8b9-7634-4312-a5e3-0c4cd75bf37f",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "k = 0 , def = 0 : 30.00 - 0.00e+00j\n",
+      "k = 0 , def = 1 : 30.00 - 0.00e+00j\n",
+      "k = 0 , def = 2 : 30.00 - 0.00e+00j\n",
+      "k = 0 , def = 3 : 30.00 - 0.00e+00j\n",
+      "k = 0 , def = 4 : 30.00 - 0.00e+00j\n",
+      "k = 1 , def = 0 : -8.00 - 0.00e+00j\n",
+      "k = 1 , def = 1 : -8.00 + 2.66e-15j\n",
+      "k = 1 , def = 2 : -8.00 - 1.68e-15j\n",
+      "k = 1 , def = 3 : -8.00 + 1.78e-15j\n",
+      "k = 1 , def = 4 : -8.00 - 2.66e-15j\n",
+      "k = 2 , def = 0 : 6.00 - 0.00e+00j\n",
+      "k = 2 , def = 1 : 6.00 + 4.44e-15j\n",
+      "k = 2 , def = 2 : 6.00 - 2.48e-15j\n",
+      "k = 2 , def = 3 : 6.00 - 4.44e-16j\n",
+      "k = 2 , def = 4 : 6.00 - 4.44e-15j\n",
+      "k = 3 , def = 0 : -2.00 - 0.00e+00j\n",
+      "k = 3 , def = 1 : -2.00 + 7.53e-15j\n",
+      "k = 3 , def = 2 : -2.00 - 4.59e-15j\n",
+      "k = 3 , def = 3 : -2.00 + 6.61e-15j\n",
+      "k = 3 , def = 4 : -2.00 - 7.53e-15j\n",
+      "k = 4 , def = 0 : 6.00 - 0.00e+00j\n",
+      "k = 4 , def = 1 : 6.00 + 7.99e-15j\n",
+      "k = 4 , def = 2 : 6.00 - 4.07e-15j\n",
+      "k = 4 , def = 3 : 6.00 - 8.88e-16j\n",
+      "k = 4 , def = 4 : 6.00 - 7.99e-15j\n",
+      "k = 5 , def = 0 : -8.00 - 0.00e+00j\n",
+      "k = 5 , def = 1 : -8.00 + 1.47e-14j\n",
+      "k = 5 , def = 2 : -8.00 - 9.76e-15j\n",
+      "k = 5 , def = 3 : -8.00 + 1.20e-14j\n",
+      "k = 5 , def = 4 : -8.00 - 1.47e-14j\n"
+     ]
+    }
+   ],
+   "source": [
+    "x = [4, 3, 5, 10, 5, 3]\n",
+    "\n",
+    "def format_complex(z):\n",
+    "\tpm = \" + \" if z.imag < 0 else \" - \"\n",
+    "\timag = f\"{np.abs(z.imag):.2e}\" if np.abs(z.imag) < 1e-6 else f\"{np.abs(z.imag):.3f}\"\n",
+    "\treturn f\"{z.real:.2f}\" + pm + imag + \"j\"\n",
+    "\n",
+    "for k in range(0, len(x)):\n",
+    "\tfor i,f in enumerate([FFT, DFT, DFT_, DFT__, DFT_even]):\n",
+    "\t\tprint(\"k =\", k, \", def =\", i+1, \":\", format_complex(f(k)))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "4be9a3ba-49db-4f62-824c-f32705ff9e12",
+   "metadata": {},
+   "source": [
+    "Note that if we do it with a signal that *isn't* perfectly even, the last definition doesn't equal the others."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "id": "bd707baf-e504-476a-a038-d85bc261ce64",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "k = 0 , def = 1 : 29.00 - 0.00e+00j\n",
+      "k = 0 , def = 2 : 29.00 - 0.00e+00j\n",
+      "k = 0 , def = 5 : 29.00 - 0.00e+00j\n",
+      "k = 1 , def = 1 : -8.50 + 0.866j\n",
+      "k = 1 , def = 2 : -8.50 + 0.866j\n",
+      "k = 1 , def = 5 : -8.50 - 0.866j\n",
+      "k = 2 , def = 1 : 6.50 + 0.866j\n",
+      "k = 2 , def = 2 : 6.50 + 0.866j\n",
+      "k = 2 , def = 5 : 6.50 - 0.866j\n",
+      "k = 3 , def = 1 : -1.00 - 1.11e-16j\n",
+      "k = 3 , def = 2 : -1.00 + 5.14e-15j\n",
+      "k = 3 , def = 5 : -1.00 - 5.14e-15j\n",
+      "k = 4 , def = 1 : 6.50 - 0.866j\n",
+      "k = 4 , def = 2 : 6.50 - 0.866j\n",
+      "k = 4 , def = 5 : 6.50 + 0.866j\n",
+      "k = 5 , def = 1 : -8.50 - 0.866j\n",
+      "k = 5 , def = 2 : -8.50 - 0.866j\n",
+      "k = 5 , def = 5 : -8.50 + 0.866j\n"
+     ]
+    }
+   ],
+   "source": [
+    "x = [4, 3, 5, 10, 5, 2]\n",
+    "\n",
+    "for k in range(0, len(x)):\n",
+    "\tfor i,f in zip((1, 2, 5), [FFT, DFT, DFT_even]):\n",
+    "\t\tprint(\"k =\", k, \", def =\", i, \":\", format_complex(f(k)))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "7cdc97c4-3b5a-4539-a45d-db19d79ff378",
+   "metadata": {},
+   "source": [
+    "[\"One, two, five--\"  \n",
+    "\"Three, sir!\"  \n",
+    "\"Three!\"](https://www.youtube.com/watch?v=xOrgLj9lOwk&t=129s)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "84cfeed2-f49c-4b46-a6db-c5d5cf94f8a2",
+   "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.13.1"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

specderiv/__init__.py

@@ -1,3 +1,3 @@
 from .specderiv import cheb_deriv, fourier_deriv
 
-__version__ = '0.4'
+__version__ = '0.5'