diff --git a/org/qmckl_jastrow.org b/org/qmckl_jastrow.org
index c017e8f..ed7de29 100644
--- a/org/qmckl_jastrow.org
+++ b/org/qmckl_jastrow.org
@@ -2530,7 +2530,7 @@ qmckl_exit_code qmckl_provide_factor_en_deriv_e(qmckl_context context)
   return QMCKL_SUCCESS;
 }
     #+end_src
-
+    
 *** Compute
    :PROPERTIES:
    :Name:     qmckl_compute_factor_en_deriv_e
@@ -3142,9 +3142,10 @@ assert(fabs(een_rescaled_e[0][1][5][2]-0.3424402276009091)    < 1.e-12);
 
 ** Electron-electron rescaled distances for each order and derivatives
 
-   ~een_rescaled_e~ stores the table of the rescaled distances between all
-   pairs of electrons and raised to the power \(p\) defined by ~cord_num~.
-   Here we take its derivatives required for the een jastrow.
+   ~een_rescaled_e_deriv_e~ stores the table of the derivatives of the
+   rescaled distances between all pairs of electrons and raised to the
+   power \(p\) defined by ~cord_num~.  Here we take its derivatives
+   required for the een jastrow.
    
    TODO: write formulae
 
@@ -3414,6 +3415,16 @@ for i in range(elec_num):
   for j in range(elec_num):
     elec_dist[i, j] = np.linalg.norm(elec_coord[i] - elec_coord[j])
 
+elec_dist_deriv_e = np.zeros(shape=(4,elec_num, elec_num),dtype=float)
+for j in range(elec_num):
+  for i in range(elec_num):
+    rij_inv = 1.0 / elec_dist[i, j]
+    for ii in range(3):
+      elec_dist_deriv_e[ii, i, j] = -(elec_coord[j][ii] - elec_coord[i][ii]) * rij_inv
+    elec_dist_deriv_e[3, i, j] = 2.0 * rij_inv
+  elec_dist_deriv_e[:, j, j] = 0.0
+
+
 kappa = 1.0
 
 een_rescaled_e_ij = np.zeros(shape=(elec_num * (elec_num - 1)//2, cord_num+1), dtype=float)
@@ -3441,25 +3452,49 @@ for l in range(1,cord_num+1):
       een_rescaled_e[j, i, l] = x
       k = k + 1
 
-print(" een_rescaled_e[0, 2, 1] = ",een_rescaled_e[0, 2, 1])
-print(" een_rescaled_e[0, 3, 1] = ",een_rescaled_e[0, 3, 1])
-print(" een_rescaled_e[0, 4, 1] = ",een_rescaled_e[0, 4, 1])
-print(" een_rescaled_e[1, 3, 2] = ",een_rescaled_e[1, 3, 2])
-print(" een_rescaled_e[1, 4, 2] = ",een_rescaled_e[1, 4, 2])
-print(" een_rescaled_e[1, 5, 2] = ",een_rescaled_e[1, 5, 2])
+een_rescaled_e_deriv_e = np.zeros(shape=(elec_num,4,elec_num,cord_num+1),dtype=float)
+for l in range(0,cord_num+1):
+  kappa_l = -1.0 * kappa * l
+  for j in range(0,elec_num):
+    for i in range(0,elec_num):
+      for ii in range(0,4):
+        een_rescaled_e_deriv_e[i,ii,j,l] = kappa_l * elec_dist_deriv_e[ii,i,j] 
+      een_rescaled_e_deriv_e[i,3,j,l] = een_rescaled_e_deriv_e[i,3,j,l] +        \
+      een_rescaled_e_deriv_e[i,0,j,l] * een_rescaled_e_deriv_e[i,0,j,l] +        \
+      een_rescaled_e_deriv_e[i,1,j,l] * een_rescaled_e_deriv_e[i,1,j,l] +        \
+      een_rescaled_e_deriv_e[i,2,j,l] * een_rescaled_e_deriv_e[i,2,j,l] 
+
+      for ii in range(0,4):
+        een_rescaled_e_deriv_e[i,ii,j,l] = een_rescaled_e_deriv_e[i,ii,j,l] * een_rescaled_e[i,j,l]
+
+print(" een_rescaled_e_deriv_e[1, 1, 3, 1] = ",een_rescaled_e_deriv_e[0, 0, 2, 1])
+print(" een_rescaled_e_deriv_e[1, 1, 4, 1] = ",een_rescaled_e_deriv_e[0, 0, 3, 1])
+print(" een_rescaled_e_deriv_e[1, 1, 5, 1] = ",een_rescaled_e_deriv_e[0, 0, 4, 1])
+print(" een_rescaled_e_deriv_e[2, 1, 4, 2] = ",een_rescaled_e_deriv_e[1, 0, 3, 2])
+print(" een_rescaled_e_deriv_e[2, 1, 5, 2] = ",een_rescaled_e_deriv_e[1, 0, 4, 2])
+print(" een_rescaled_e_deriv_e[2, 1, 6, 2] = ",een_rescaled_e_deriv_e[1, 0, 5, 2])
     #+end_src
 
     #+RESULTS:
-    :  een_rescaled_e[0, 2, 1] =  0.08084493981483197
-    :  een_rescaled_e[0, 3, 1] =  0.1066745707571846
-    :  een_rescaled_e[0, 4, 1] =  0.017542731694647366
-    :  een_rescaled_e[1, 3, 2] =  0.02214680362033448
-    :  een_rescaled_e[1, 4, 2] =  0.0005700154999202759
-    :  een_rescaled_e[1, 5, 2] =  0.3424402276009091
+    :  een_rescaled_e_deriv_e[1, 1, 3, 1] =  0.05991352796887283
+    :  een_rescaled_e_deriv_e[1, 1, 4, 1] =  0.011714035071545248
+    :  een_rescaled_e_deriv_e[1, 1, 5, 1] =  0.00441398875758468
+    :  een_rescaled_e_deriv_e[2, 1, 4, 2] =  0.013553180060167595
+    :  een_rescaled_e_deriv_e[2, 1, 5, 2] =  0.00041342909359870457
+    :  een_rescaled_e_deriv_e[2, 1, 6, 2] =  0.5880599146214673
 
      #+begin_src c :tangle (eval c_test)
 //assert(qmckl_electron_provided(context));
+double een_rescaled_e_deriv_e[walk_num][elec_num][4][elec_num][(cord_num + 1)];
+rc = qmckl_get_jastrow_een_rescaled_e_deriv_e(context, &(een_rescaled_e_deriv_e[0][0][0][0][0]));
 
+// value of (0,0,0,2,1)
+assert(fabs(een_rescaled_e_deriv_e[0][0][0][2][1] + 0.05991352796887283   ) < 1.e-12);
+assert(fabs(een_rescaled_e_deriv_e[0][0][0][3][1] + 0.011714035071545248  ) < 1.e-12);
+assert(fabs(een_rescaled_e_deriv_e[0][0][0][4][1] + 0.00441398875758468   ) < 1.e-12);
+assert(fabs(een_rescaled_e_deriv_e[0][1][0][3][2] + 0.013553180060167595  ) < 1.e-12);
+assert(fabs(een_rescaled_e_deriv_e[0][1][0][4][2] + 0.00041342909359870457) < 1.e-12);
+assert(fabs(een_rescaled_e_deriv_e[0][1][0][5][2] + 0.5880599146214673    ) < 1.e-12);
      #+end_src
 
 ** Electron-nucleus rescaled distances for each order
@@ -4092,15 +4127,15 @@ print(" een_rescaled_n[1, 5, 2] = ",een_rescaled_n[1, 5, 2])
      #+begin_src c :tangle (eval c_test)
 assert(qmckl_electron_provided(context));
 
-//double een_rescaled_n_deriv_e[walk_num][elec_num][4][nucl_num][(cord_num + 1)];
-//rc = qmckl_get_jastrow_een_rescaled_n_deriv_e(context, &(een_rescaled_n_deriv_e[0][0][0][0][0]));
+double een_rescaled_n_deriv_e[walk_num][elec_num][4][nucl_num][(cord_num + 1)];
+rc = qmckl_get_jastrow_een_rescaled_n_deriv_e(context, &(een_rescaled_n_deriv_e[0][0][0][0][0]));
 
-//printf(" 0 2 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][2][0][0][1]);
-//printf(" 0 3 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][3][0][0][1]);
-//printf(" 0 4 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][4][0][0][1]);
-//printf(" 0 3 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][3][0][1][2]);
-//printf(" 0 4 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][4][0][1][2]);
-//printf(" 0 5 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][5][0][1][2]);
+printf(" 0 2 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][2][0][0][1]);
+printf(" 0 3 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][3][0][0][1]);
+printf(" 0 4 0 1 %10.16f\n",een_rescaled_n_deriv_e[0][4][0][0][1]);
+printf(" 0 3 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][3][0][1][2]);
+printf(" 0 4 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][4][0][1][2]);
+printf(" 0 5 1 2 %10.16f\n",een_rescaled_n_deriv_e[0][5][0][1][2]);
 // value of (0,2,1)
 //assert(fabs(een_rescaled_n_deriv_e[0][2][0][1]-0.10612983920006765)  < 1.e-12);
 //assert(fabs(een_rescaled_n_deriv_e[0][3][0][1]-0.135652809635553)    < 1.e-12);