From 1ed479e965158551f0b13faf9d3d31330c5184f7 Mon Sep 17 00:00:00 2001
From: mcarfagno <m.carfagno@outlook.com>
Date: Mon, 14 Dec 2020 14:09:49 +0000
Subject: [PATCH] slight rework of numerical jacobian functions

---
 notebooks/numericalJacobian.ipynb | 52 +++++++++++++++----------------
 1 file changed, 25 insertions(+), 27 deletions(-)

diff --git a/notebooks/numericalJacobian.ipynb b/notebooks/numericalJacobian.ipynb
index dd60b5d..cd28bbf 100644
--- a/notebooks/numericalJacobian.ipynb
+++ b/notebooks/numericalJacobian.ipynb
@@ -15,7 +15,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 31,
+   "execution_count": 19,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -69,7 +69,7 @@
     "        theta + v*np.arctan(delta)/L*dt\n",
     "        ])\n",
     "\n",
-    "def J(f,x,u,epsilon=1e-4):\n",
+    "def Jacobians(f,x,u,epsilon=1e-4):\n",
     "    \"\"\"\n",
     "    :param f:\n",
     "    :param x:\n",
@@ -77,31 +77,36 @@
     "    \"\"\"\n",
     "    \n",
     "    jac_x = np.zeros((4,4))\n",
-    "    \n",
-    "    for i in range(x.shape[0]):\n",
-    "        perturb = np.zeros(x.shape[0])\n",
-    "        perturb[i] = epsilon\n",
-    "        \n",
-    "        #each loumn of the jac is given by perturbing a variable\n",
-    "        jac_x[:,i]= (f(x+perturb, u)-f(x-perturb, u))/2*epsilon\n",
-    "    \n",
-    "    \n",
     "    jac_u = np.zeros((4,2))\n",
     "    \n",
+    "    perturb_x = np.eye(4)*epsilon\n",
+    "    perturb_u = np.eye(2)*epsilon\n",
+    "    \n",
+    "    #each row is state vector where one variable has been perturbed\n",
+    "    x_perturbed_plus = np.tile(x,(4,1))+perturb_x\n",
+    "    x_perturbed_minus = np.tile(x,(4,1))-perturb_x\n",
+    "    \n",
+    "    #each row is state vector where one variable has been perturbed\n",
+    "    u_perturbed_plus = np.tile(u,(2,1))+perturb_u\n",
+    "    u_perturbed_minus = np.tile(u,(2,1))-perturb_u\n",
+    "    \n",
+    "    for i in range(x.shape[0]):\n",
+    "        \n",
+    "        #each coloumn of the jac is given by perturbing a variable\n",
+    "        jac_x[:,i]= (f(x+perturb_x[i,:], u)-f(x-perturb_x[i,:], u))/2*epsilon\n",
+    "        \n",
     "    for i in range(u.shape[0]):\n",
-    "        perturb = np.zeros(u.shape[0])\n",
-    "        perturb[i] = epsilon\n",
-    "        \n",
-    "        #each loumn of the jac is given by perturbing a variable\n",
-    "        jac_u[:,i]= (f(x, u+perturb)-f(x, u-perturb))/2*epsilon\n",
     "        \n",
+    "        #each coloumn of the jac is given by perturbing a variable\n",
+    "        jac_u[:,i]= (f(x, u+perturb_u[i,:])-f(x, u-perturb_u[i,:]))/2*epsilon\n",
+    "\n",
     "    return jac_x, jac_u\n",
     "    "
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 32,
+   "execution_count": 20,
    "metadata": {},
    "outputs": [
     {
@@ -117,7 +122,7 @@
        "        [0.0000000e+00, 3.2051282e-09]]))"
       ]
      },
-     "execution_count": 32,
+     "execution_count": 20,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -127,15 +132,8 @@
     "x=np.array([0,0,1,0])\n",
     "u=np.array([1,0.2])\n",
     "\n",
-    "J(f,x,u)"
+    "Jacobians(f,x,u)"
    ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": []
   }
  ],
  "metadata": {
@@ -154,7 +152,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.7.6"
+   "version": "3.8.5"
   }
  },
  "nbformat": 4,