Add some maths

experiments/ELM.ipynb

@@ -2,7 +2,7 @@
  "cells": [
   {
    "cell_type": "code",
-   "execution_count": 1,
+   "execution_count": 6,
    "id": "e465f900-0a94-46d4-b3c5-23cc0129557b",
    "metadata": {},
    "outputs": [],
@@ -16,7 +16,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 72,
+   "execution_count": 7,
    "id": "b39df788-2971-4974-9ced-1518a47181e8",
    "metadata": {},
    "outputs": [],
@@ -33,7 +33,7 @@
     "    # c\n",
     "    bond_bending_constant: float\n",
     "    \n",
-    "    def make(n_rows: int, n_cols: int) -> ELM:\n",
+    "    def make(n_rows: int, n_cols: int):\n",
     "        return ELM(\n",
     "            n_rows,\n",
     "            n_cols,\n",
@@ -118,7 +118,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 77,
+   "execution_count": 8,
    "id": "b4ed52bd-02aa-411a-9af4-be9bb28227bc",
    "metadata": {},
    "outputs": [
@@ -207,7 +207,38 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 82,
+   "execution_count": null,
+   "id": "0f1f9dfd-8266-4874-9ba5-83a8a784d0e9",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import multiprocessing\n",
+    "\n",
+    "def velocity_verlet_step(elm: ELM, delta_t: float) -> ELM:\n",
+    "    new_displacement = elm.displacement + elm.velocity * delta_t + 0.5 * elm.force * delta_t * delta_t\n",
+    "    \n",
+    "    new_force = numpy.zeros((elm.n_rows, elm.n_cols, 2))\n",
+    "\n",
+    "    def job_force(row):\n",
+    "        for col in range(elm.n_cols):\n",
+    "            new_force[row, col, :] = calculate_force(\n",
+    "                elm,\n",
+    "                new_displacement,\n",
+    "                row,\n",
+    "                col\n",
+    "            )\n",
+    "\n",
+    "    with multiprocessing.Pool(4) as pool:\n",
+    "        pool.map(new_force, range(elm.n_rows))\n",
+    "    \n",
+    "    new_velocity = elm.velocity + 0.5 * (elm.force + new_force) * delta_t\n",
+    "\n",
+    "    return ELM(elm.n_rows, elm.n_cols, new_displacement, new_velocity, new_force)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 9,
    "id": "ff7adc30-191e-4b7c-b146-2a533463d6ea",
    "metadata": {},
    "outputs": [],
@@ -232,7 +263,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 85,
+   "execution_count": 10,
    "id": "bea964bf-6e94-4b72-a458-50c2964bbba1",
    "metadata": {},
    "outputs": [],
@@ -244,14 +275,14 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 86,
+   "execution_count": 12,
    "id": "421c0fd2-2dba-4fd5-89ff-0968227b08b9",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "0e9725fa01fe4ff6a520690877b4f5eb",
+       "model_id": "7d7f0752567340c08c97df22788fc638",
        "version_major": 2,
        "version_minor": 0
       },
@@ -268,7 +299,7 @@
        "<function __main__.do_plot(i: int) -> None>"
       ]
      },
-     "execution_count": 86,
+     "execution_count": 12,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -297,7 +328,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "38e341e1-de34-41e7-b8e7-1322d27db18b",
+   "id": "96ab8cf1-8739-4faa-b6b5-70aa52a15a2a",
    "metadata": {},
    "outputs": [],
    "source": []
@@ -319,7 +350,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.10.11"
+   "version": "3.9.13"
   }
  },
  "nbformat": 4,

experiments/ELM_Debug.ipynb

@@ -64,7 +64,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 18,
+   "execution_count": 4,
    "id": "82bb4a02-1ecf-4c84-8b9d-b753dee75337",
    "metadata": {},
    "outputs": [
@@ -78,10 +78,10 @@
     {
      "data": {
       "text/plain": [
-       "array([-0.14798801, -0.13453456])"
+       "0"
       ]
      },
-     "execution_count": 18,
+     "execution_count": 4,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -99,6 +99,7 @@
     "\n",
     "    print(f\"{j=} {K_ij=} {u_ij=} {x_ij=}\")\n",
     "\n",
+    "    return 0\n",
     "    return -K_ij * np.dot(u_ij, x_ij) * r_ij / np.linalg.norm(r_ij)\n",
     "\n",
     "def calculate_force_bond(\n",
@@ -129,14 +130,14 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 19,
+   "execution_count": 8,
    "id": "f7a019c2-8d20-4e84-9e0c-ef788781a07a",
    "metadata": {},
    "outputs": [
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "78145432b809422bb171f675b00d05a6",
+       "model_id": "10515f6b18554cd1bd86dd56891eed2b",
        "version_major": 2,
        "version_minor": 0
       },
@@ -153,7 +154,7 @@
        "<function __main__.do_plot(dp_x: float, dp_y: float) -> None>"
       ]
      },
-     "execution_count": 19,
+     "execution_count": 8,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -191,9 +192,9 @@
     "        [pi[1]] * len(elastic_ys),\n",
     "        elastic_xs,\n",
     "        elastic_ys,\n",
-    "        scale = 5,\n",
     "        angles = \"xy\",\n",
     "        pivot = \"tip\",\n",
+    "        scale = 5,\n",
     "        width = 0.005\n",
     "    )\n",
     "    plt.quiver(\n",
@@ -247,7 +248,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.10.11"
+   "version": "3.9.13"
   }
  },
  "nbformat": 4,